ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/REPOS_ERICCA/magic/lib/CAD4FE/src/CAD4FE_LocalEdgeCriteria.cpp
Revision: 569
Committed: Thu Oct 16 14:36:31 2014 UTC (10 years, 6 months ago) by foucault
File size: 16712 byte(s)
Error occurred while calculating annotation data.
Log Message:
Mise à jour pour CAD4FE (Gilles) : operation 1 (tentative)

File Contents

# Content
1 //---------------------------------------------------------------------------
2 #include <fstream>
3 #include <sstream>
4
5 //---------------------------------------------------------------------------
6 // include MAGIC Headers
7 #include "gestionversion.h"
8 #include <mg_geometrie.h>
9 #include <mg_maillage.h>
10 #include <mg_maillage_outils.h>
11 #include <ot_mathematique.h>
12
13 #include "CAD4FE_MCEdge.h"
14 #include "CAD4FE_MCAA.h"
15 #include "CAD4FE_MCBody.h"
16 #include "CAD4FE_mg_utils.h"
17 #include "CAD4FE_Intersection_Plane_MG_MAILLAGE.h"
18 #include "ot_algorithme_geometrique.h"
19 #include "CAD4FE_Criteria.h"
20 #include "CAD4FE_ColorMap.h"
21
22 #pragma hdrstop
23
24 #include "CAD4FE_LocalEdgeCriteria.h"
25
26 //---------------------------------------------------------------------------
27
28 #pragma package(smart_init)
29
30 using namespace CAD4FE;
31
32 //------------------------------------------------------------------------------
33 // Edge Local Properties
34
35 LocalEdgeCriteria::LocalEdgeCriteria(MG_SEGMENT * __seg, MCAA * __mcaa, MG_MAILLAGE * __mesh, MG_SEGMENT * __startSeg)
36 : _seg(__seg), _mcaa(__mcaa), _mesh(__mesh), _startSeg(__startSeg), _normalOffset(NULL)
37 {
38 double x[2][3];
39 MG_UTILS::MG_NOEUD_GET_XYZ(_seg->get_noeud1(), x[0]);
40 MG_UTILS::MG_NOEUD_GET_XYZ(_seg->get_noeud2(), x[1]);
41 // _deletionScoreOppositeEdge = 0;
42
43 /* GF Debug int nb_triangles_adj_start_segment = __startSeg->get_lien_triangle()->get_nb();
44 printf("start seg %ld : %d triangles adjacents\n", __startSeg->get_id(), nb_triangles_adj_start_segment);*/
45
46 for (unsigned i=0; i<3; i++)
47 {
48 _P[i] = .5*(x[0][i]+x[1][i]);
49 _N[i] = x[1][i]-_P[i];
50 }
51
52 _meshSize = 0;
53 _normalOffset = 0;
54
55 _mcEdge = NULL;
56
57 Update();
58 if (_mcEdge) time = _mcEdge->time;
59 }
60 //------------------------------------------------------------------------------
61
62 LocalEdgeCriteria::~LocalEdgeCriteria()
63 {
64 delete _normalOffset;
65 }
66 //------------------------------------------------------------------------------
67
68 void
69 LocalEdgeCriteria::Update(bool __reconstructNormalOffset)
70 {
71 _meshSize = _mcaa->GetSize(_P);
72
73 double targetLength = .86*_meshSize;
74 if ( __reconstructNormalOffset
75 ||!_normalOffset )
76 // inutile? || fabs ( 1 - 2*targetLength/_normalOffset->GetTargetLength()) > .1)
77 {
78 if (_normalOffset) delete _normalOffset;
79 _normalOffset = new Intersection_Plane_MG_MAILLAGE(_P, _N, _mesh, targetLength, _startSeg, NULL);
80 }
81
82 InitSetOfTouchingEdges();
83
84 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
85 int nb_pts = pl.size();
86 int o = _normalOffset->GetOriginIndex();
87
88 // Get the mc edge attached to the start segment
89 {
90 MG_SEGMENT * seg = pl[o].E;
91 MG_ELEMENT_TOPOLOGIQUE * topo = seg->get_lien_topologie();
92 if ( topo && topo->get_dimension() == 1)
93 {
94 MCEdge * edge = (MCEdge*) topo;
95 _mcEdge = edge;
96 time = _mcEdge->time;
97 }
98 else
99 {
100 _mcEdge = (MCEdge*)_seg->get_lien_topologie();
101 }
102 }
103
104 for (unsigned i=0; i<2; i++)
105 {
106 _lengthToOppositeEdges[i] = 0;
107 _oppositeEdgesPolylineIndices[i] = -1;
108 _oppositeEdges[i] = NULL;
109 }
110
111 for (int j = 0; j < 2; j++)
112 {
113 for (int i = 0; i+1 < nb_pts; i++)
114 {
115 int k[2];
116 if (j == 0)
117 {
118 k[0] = o - i;
119 k[1] = o - i - 1;
120 if (k[1] < 0 || k[1] >= nb_pts)
121 break;
122 }
123 else if (j == 1)
124 {
125 k[0] = o + i;
126 k[1] = o + i + 1;
127 if (k[1] < 0 || k[1] >= nb_pts)
128 break;
129 }
130
131 MG_SEGMENT * Seg [2];
132 MG_NOEUD * No [2];
133 OT_VECTEUR_3D xyz[2];
134 double dl;
135
136 for (unsigned l = 0; l < 2; l++)
137 {
138 xyz[l] = pl[k[l]].X;
139 Seg[l] = pl[k[l]].E;
140 No[l] = pl[k[l]].V;
141 }
142
143 MG_ELEMENT_TOPOLOGIQUE * topo;
144 topo = Seg[1]->get_lien_topologie();
145
146 if (topo && topo->get_dimension() == 1 )
147 {
148 _oppositeEdges[j] = (MCEdge *) topo;
149 _oppositeEdgesPolylineIndices[j] = k[1];
150 }
151 else
152 {
153 _oppositeEdges[j] = NULL;
154 _oppositeEdgesPolylineIndices[j] = -1;
155 }
156
157 dl = (xyz[0] - xyz[1]).get_longueur();
158 _lengthToOppositeEdges[j] += dl;
159
160 if ( _oppositeEdges[j] )
161 break;
162 }
163 }
164
165 _point = pl[o].X;
166 double * segStart = pl[0].X, * segEnd = pl[nb_pts-1].X;
167 for (unsigned i=0;i<3;i++)
168 {
169 _segment [0][i] = segStart[i];
170 _segment [1][i] = segEnd[i];
171 }
172
173 /** Compute the local edge properties */
174 // Face width
175 _faceWidth = std::min( _lengthToOppositeEdges[0], _lengthToOppositeEdges[1] );
176 // Epsilon (distance between the mesh and the geometry
177 _epsilon = OT_ALGORITHME_GEOMETRIQUE::Dist3D_Point_Segment( _segment [0], _segment [1], _point); // Angle between the 2 edge segments
178 if ( _segment [0] == _point || _segment [1] == _point )
179 _deviationAngle = 0;
180 else
181 _deviationAngle = OT_ALGORITHME_GEOMETRIQUE::Angle3D_Segment_Segment( _segment [0], _point, _segment [1] );
182 }
183 //------------------------------------------------------------------------------
184
185 double LocalEdgeCriteria::DeletionScore_FaceWidth() const
186 {
187 double criterionFaceWidth;
188 double limitFaceWidth = _meshSize / _mcaa->GetMaxOverdensity();
189
190 if ( _oppositeEdges[0] || _oppositeEdges[1] )
191 {
192 if (_faceWidth < limitFaceWidth)
193 criterionFaceWidth = 1 - _faceWidth / limitFaceWidth;
194 else
195 criterionFaceWidth = .05 * (1 - _faceWidth / limitFaceWidth);
196 }
197 else
198 criterionFaceWidth = -.05;
199
200 return criterionFaceWidth;
201 }
202 //------------------------------------------------------------------------------
203
204 double LocalEdgeCriteria::DeletionScore_DeviationAngle() const
205 {
206 double criterionAngle;
207 double limitAngle = _mcaa->GetLimitAngle();
208
209 if (_deviationAngle < limitAngle)
210 criterionAngle = 1 - _deviationAngle / limitAngle;
211 else
212 criterionAngle = .2*(1 - _deviationAngle / limitAngle);
213
214 return criterionAngle;
215 }
216 //------------------------------------------------------------------------------
217
218 double LocalEdgeCriteria::DeletionScore() const
219 {
220 int i, NB_CRITERIA = 0;
221 double score = 0;
222 double val[10];
223
224 if (_mcaa->GetMCBody()->G21()->GetArc(_mcEdge->get_id())->Rank() != 2)
225 return 0;
226
227 if (_mcEdge && _mcEdge->get_nb_ccf())
228 return 0;
229
230 val[NB_CRITERIA++] = DeletionScore_FaceWidth();
231 val[NB_CRITERIA++] = DeletionScore_DeviationAngle();
232 // val[NB_CRITERIA++] = DeletionScore_OppositeEdge();
233
234 for (i=0;i<NB_CRITERIA;i++)
235 {
236 if (score < val[i])
237 {
238 score = val[i];
239 }
240 }
241
242 return score;
243 }
244 //------------------------------------------------------------------------------
245
246 void LocalEdgeCriteria::GetPolylineVertex(int i, float * __vec3f) const
247 {
248 const double * vec3d = _normalOffset->GetPolyline()[i].X;
249 for (int i=0; i<3; i++)
250 __vec3f[i] = vec3d[i];
251 }
252 //------------------------------------------------------------------------------
253
254 unsigned LocalEdgeCriteria::GetPolylineVerticesCount() const
255 {
256 return _normalOffset->GetPolyline().size();
257 }
258 //------------------------------------------------------------------------------
259
260 std::string
261 LocalEdgeCriteria::InventorText()
262 {
263 std::ostringstream out;
264
265 unsigned char rgb[3];
266 ColorMap::jetColorMap(rgb, 1-DeletionScore(), 0, 1);
267 _normalOffset->SetColor(rgb);
268 out << _normalOffset->InventorText();
269
270 for (unsigned j=0; j<2; j++)
271 if (_oppositeEdgesPolylineIndices[j] != -1) {
272 int i = _oppositeEdgesPolylineIndices[j];
273 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
274 out << "\nSeparator { #sep1 \n";
275 out << "\n Coordinate3 {\n point [ \n";
276 out << pl[i].X[0] << " " <<pl[i].X[1] << " " << pl[i].X[2] << " \n";
277 out << "\n]\n}\n";
278 out << "PolygonOffset { \n";
279 out << "styles POINTS \n";
280 out << "factor 4.0 \n";
281 out << "units 1.0 \n";
282 out << "} \n";
283 out << "DrawStyle {\npointSize 4\n}\n";
284 out << "BaseColor { \n rgb 1.0 0.0 0.0\n }\n";
285 out << "PointSet {\nstartIndex "<<0<<"\nnumPoints "<<1<<"\n}\n";
286 out << "} # end Sep1 \n";
287 }
288
289 return out.str();
290 }
291 //------------------------------------------------------------------------------
292
293 MCEdge * LocalEdgeCriteria::GetEdge()
294 {
295 return _mcEdge;
296 }
297 //------------------------------------------------------------------------------
298
299 MG_SEGMENT * LocalEdgeCriteria::GetSegment()
300 {
301 return _seg;
302 }
303 //------------------------------------------------------------------------------
304
305 double LocalEdgeCriteria::GetFaceWidth()
306 {
307 return _faceWidth;
308 }
309 //------------------------------------------------------------------------------
310
311 double LocalEdgeCriteria::GetEpsilon()
312 {
313 return _epsilon;
314 }
315 //------------------------------------------------------------------------------
316
317 double LocalEdgeCriteria::GetDeviationAngle()
318 {
319 return _deviationAngle;
320 }
321 //------------------------------------------------------------------------------
322 bool LocalEdgeCriteria::IsTouchingEdge(MCEdge * __mcEdge)
323 {
324 return _setTouchingEdge.find(__mcEdge) != _setTouchingEdge.end();
325 }
326 //------------------------------------------------------------------------------
327 bool LocalEdgeCriteria::IsTouchingSegment(MG_SEGMENT * __segment)
328 {
329 if (__segment == NULL)
330 return false;
331
332 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
333 unsigned nb_pts = pl.size();
334
335 // Get the mc edge attached to the start segment
336 for (unsigned i = 0; i < nb_pts; i++)
337 {
338 MG_SEGMENT * seg = pl[i].E;
339
340 if (seg == __segment)
341 return true;
342 }
343
344 return false;
345 }
346 //------------------------------------------------------------------------------
347 bool LocalEdgeCriteria::IsStartSegment(MG_SEGMENT * __segment)
348 {
349 if (__segment == NULL)
350 return false;
351
352 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
353 unsigned int o = _normalOffset->GetOriginIndex();
354
355 return (pl[o].E == __segment);
356 }
357 //------------------------------------------------------------------------------
358
359 MCEdge * LocalEdgeCriteria::GetTouchingEdge(int __index)
360 {
361 if (__index >= 0 && __index < 2)
362 return _oppositeEdges[__index];
363 else
364 {
365 printf("Error GetTouchingEdge: __index = %d\n", __index);
366 return NULL;
367 }
368 }
369 //------------------------------------------------------------------------------
370
371 double LocalEdgeCriteria::GetLengthToTouchingEdge(int __index)
372 {
373 if (__index >= 0 && __index < 2)
374 return _lengthToOppositeEdges[__index];
375 else
376 {
377 printf("Error GetLengthToTouchingEdge: __index = %d\n", __index);
378 return -1;
379 }
380 }
381 //------------------------------------------------------------------------------
382 MCEdge * LocalEdgeCriteria::GetClosestTouchingEdge()
383 {
384 return (_lengthToOppositeEdges[0] < _lengthToOppositeEdges[1])
385 ? _oppositeEdges[0]
386 : _oppositeEdges[1];
387 }
388 //------------------------------------------------------------------------------
389 double * LocalEdgeCriteria::GetClosestTouchingEdgePoint()
390 {
391 int i = (_lengthToOppositeEdges[0] < _lengthToOppositeEdges[1])
392 ? 0
393 : 1;
394 return _normalOffset->GetPolyline()[ _oppositeEdgesPolylineIndices [i] ].X;
395 }
396
397 double * LocalEdgeCriteria::GetTouchingEdgePoint(int __index)
398 {
399 if (__index >= 0 && __index < 2)
400 return _normalOffset->GetPolyline()[ _oppositeEdgesPolylineIndices [__index] ].X;
401 else
402 {
403 printf("Error GetLengthToTouchingEdge: __index = %d\n", __index);
404 return 0;
405 }
406 }
407 //------------------------------------------------------------------------------
408 void LocalEdgeCriteria::CheckMCTess()
409 {
410 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
411 int nb_pts = pl.size();
412 int o = _normalOffset->GetOriginIndex();
413
414 for (int j = 0; j < 2; j++)
415 {
416 for (int i = 0; i+1 < nb_pts; i++)
417 {
418 int k[2];
419 if (j == 0)
420 {
421 k[0] = o - i;
422 k[1] = o - i - 1;
423 if (k[1] < 0 || k[1] >= nb_pts)
424 break;
425 }
426 else if (j == 1)
427 {
428 k[0] = o + i;
429 k[1] = o + i + 1;
430 if (k[1] < 0 || k[1] >= nb_pts)
431 break;
432 }
433
434 MG_SEGMENT * Seg [2];
435 MG_NOEUD * No [2];
436 OT_VECTEUR_3D xyz[2];
437 double dl;
438
439 for (unsigned l = 0; l < 2; l++)
440 {
441 xyz[l] = pl[k[l]].X;
442 Seg[l] = pl[k[l]].E;
443 No[l] = pl[k[l]].V;
444 }
445
446 MG_ELEMENT_TOPOLOGIQUE * topo;
447 topo = Seg[1]->get_lien_topologie();
448
449 if ( ! _mcaa->GetMCTess()->contient(Seg[1]))
450 printf("The MC Tessellation does not contains a segment crossed by LEC's polyline !\n");
451 _mcaa->CheckIfTopoExists(topo);
452 if (topo->get_dimension() < 2)
453 break;
454 }
455 }
456 }
457
458 // void LocalEdgeCriteria::Set_DeletionScore_OppositeEdge(double __value)
459 // {
460 // _deletionScoreOppositeEdge = __value;
461 // }
462
463 // double LocalEdgeCriteria::DeletionScore_OppositeEdge() const
464 // {
465 // if (_deletionScoreOppositeEdge>0)
466 // printf("ha!\n");
467 // return _deletionScoreOppositeEdge;
468 // }
469
470 double LocalEdgeCriteria::GetMeshSize()
471 {
472 return _meshSize;
473 }
474
475 void LocalEdgeCriteria::InitSetOfTouchingEdges()
476 {
477 std::vector <Intersection_Plane_MG_MAILLAGE::Vertex> & pl = _normalOffset->GetPolyline();
478 unsigned nb_pts = pl.size();
479
480 _setTouchingEdge.clear();
481 for (unsigned i = 0; i < nb_pts; i++)
482 {
483 MG_SEGMENT * seg = pl[i].E;
484 MG_ELEMENT_TOPOLOGIQUE * topo = seg->get_lien_topologie();
485 if (topo && topo->get_dimension() == 1)
486 _setTouchingEdge.insert((MCEdge*)topo);
487 }
488 }