Branch data Line data Source code
1 : : //- Class: FacetEvalTool
2 : : //- Description: The FacetEvalTool is a general purpose tool that uses a set
3 : : //- of facets to determine certain geometry computations. This
4 : : //- class can be used to define a surface or do fast bounding box
5 : : //- calculations, or determine point inside or outside.
6 : : //- Assumptions: It is assumed/required that the facets be continuous through
7 : : //- the set. In other words they must shared nodes and be conected
8 : : //- through each other. The algorithms assume, for efficient searching,
9 : : //- that from one facet, you can get to any other facet in the set.
10 : : //- Owner: Steven J. Owen
11 : : //- Checked by:
12 : :
13 : : #ifndef SMOOTH_FACET_EVAL_TOOL_HPP
14 : : #define SMOOTH_FACET_EVAL_TOOL_HPP
15 : :
16 : : #include "DLIList.hpp"
17 : : #include "AbstractTree.hpp"
18 : : #include "CubitFacet.hpp" //For inline function.
19 : :
20 : : template <class Y> class KDDTree;
21 : :
22 : : namespace FacetEvalToolUtils {
23 : :
24 : 10758 : template <typename R> inline R determ3(R p1, R q1, R p2, R q2, R p3, R q3) {
25 : 10758 : return (q3) * ((p2) - (p1)) + (q2) * ((p1) - (p3)) + (q1) * ((p3) - (p2));
26 : : }
27 : 4818 : template <typename R> inline R sqr(R a) { return a * a; }
28 : 0 : template <typename R> inline R cube(R a) { return sqr(a) * a; }
29 : 0 : template <typename R> inline R quart(R a) { return sqr(a) * sqr(a); }
30 : : template <typename R> inline R blend(R x) {
31 : : return -2.0 * (x) * (x) * (x) + 3.0 * (x) * (x);
32 : : }
33 : :
34 : : } // namespace FacetEvalToolUtils
35 : :
36 : : class CubitBox;
37 : : class CubitFacet;
38 : : class CubitPoint;
39 : : class CubitVector;
40 : : class CubitTransformMatrix;
41 : :
42 : : class FacetEvalTool
43 : : {
44 : : private:
45 : : static double timeGridSearch;
46 : : static double timeFacetProject;
47 : : static int numEvals;
48 : :
49 : : int toolID;
50 : : int interpOrder;
51 : : //- interpolation order 0=linear, 1=gradient, 2=quadratic, 4=quartic Bezier
52 : : DLIList<CubitFacet*> myFacetList;
53 : : //- Facets that are used for evaluation.
54 : : DLIList<CubitFacet*> markedFacetList;
55 : : //- Facets that are used for evaluation.
56 : : DLIList<CubitPoint*> myPointList;
57 : : //- point list
58 : : DLIList<CubitFacetEdge*> myEdgeList;
59 : : DLIList<DLIList<CubitFacetEdge*>*> myLoopList;
60 : : //- edges only used for bezier interp
61 : : AbstractTree <CubitFacet*> *aTree;
62 : : //- used for surfaces with many facets
63 : : CubitBox *myBBox;
64 : : //- Bounding box of facets
65 : : double myArea;
66 : : //- area of the facet surface
67 : : double compareTol;
68 : : //- for comparing bounding box
69 : : CubitFacet *lastFacet;
70 : : //- The last facet evaluated
71 : : int isFlat;
72 : : //- The surface represented by the facets is flat
73 : : CubitBoolean isParameterized;
74 : : //- Surface contains a parameterization (uv values initialized)
75 : : double minDot;
76 : : //- computed from feature angle
77 : : int output_id;
78 : :
79 : : bool have_data_to_calculate_bbox(void);
80 : :
81 : : void set_up_grid_search(double geom_factor);
82 : : //!@{
83 : : //! \brief get the closest facets using an expanding tolerance
84 : : //!
85 : : //! \param this_point The point near which we want to find facets
86 : : //! \param facet_list The list of facets near the point
87 : : //! \param tol_used The highest tolerance at which the search was performed
88 : : //!
89 : : void facets_from_search_grid( CubitVector &this_point,
90 : : DLIList<CubitFacet *> &facet_list,
91 : : double &tol_used );
92 : : //!@{
93 : : //! \brief get the closest facets using a given tolerance
94 : : //!
95 : : //! \param this_point The point near which we want to find facets
96 : : //! \param compare_tol The tolerance at which the facets were found
97 : : //! \param facet_list The list of facets near the point
98 : : //!
99 : : void facets_from_search_grid( CubitVector &this_point,
100 : : double compare_tol,
101 : : DLIList<CubitFacet *> &facet_list );
102 : : //- set up a grid search data structure if we have lots of facets
103 : :
104 : : CubitStatus get_points_from_facets(DLIList<CubitFacet*> &facet_list,
105 : : DLIList<CubitPoint*> &point_list );
106 : : //- populates the point_list with points contained by the given facets.
107 : : CubitStatus get_edges_from_facets(DLIList<CubitFacet*> &facet_list,
108 : : DLIList<CubitFacetEdge*> &edge_list );
109 : : //- populates the edge_list with edges contained by the given facets
110 : :
111 : : //- populates the loop_list with edges contained by the given edge list
112 : : void destroy_facets();
113 : : //- Destroys the facets, and points.
114 : :
115 : : void debug_draw_edges(int color = -1 );
116 : : void debug_draw_edge(CubitFacetEdge *edge, int color = -1 );
117 : : void debug_draw_vec( CubitVector & vec, int color = -1 );
118 : : void debug_draw_facet_normals( int color = -1 );
119 : : void debug_draw_point_normals( int color = -1 );
120 : : void debug_draw_bezier_edges( int color = -1 );
121 : : void debug_draw_bezier_facets( int color = -1 );
122 : : void debug_draw_bezier_facet( CubitFacet *facet, int color = -1 );
123 : : void debug_draw_line( CubitVector &begin, CubitVector &end, int color = -1 );
124 : : void debug_draw_eval_bezier_facet( CubitFacet *facet );
125 : : void debug_draw_location( CubitVector &location, int color = -1 );
126 : : void write_loops();
127 : : //- drawing functions for debugging.
128 : :
129 : : CubitStatus project_to_facets(CubitVector &this_point,
130 : : CubitBoolean trim,
131 : : CubitBoolean *outside,
132 : : CubitVector *closest_point_ptr,
133 : : CubitVector *normal_ptr);
134 : :
135 : : CubitStatus init_gradient();
136 : : //- initiaize the gradients (tangent planes at the points)
137 : :
138 : : CubitStatus init_quadrics();
139 : : //- initialize quadric coefficients at the points
140 : :
141 : : CubitStatus eval_quadratic( CubitFacet *facet,
142 : : int pt_idx,
143 : : CubitVector &eval_pt,
144 : : CubitVector &qpoint,
145 : : CubitVector &qnorm );
146 : : //- evaluate on facet based on quadratic approximation
147 : :
148 : : void on_circle( CubitVector &ptA,
149 : : CubitVector &normA,
150 : : CubitVector &ptB,
151 : : CubitVector &normB,
152 : : CubitVector &eval_pt,
153 : : CubitVector &pt_on_circle,
154 : : CubitVector &normal );
155 : : //- compute a projection of a point onto a circle. The circle
156 : : //- is defined by two points and two normals. Return the
157 : : //- projected point and its normal
158 : :
159 : : void eval_quadric( CubitFacet *facet,
160 : : int pt_index,
161 : : CubitVector &eval_pt,
162 : : CubitVector &qpt );
163 : : //- evaluate on facet based on quadratic approximation
164 : :
165 : : CubitStatus get_close_points( CubitPoint *point,
166 : : CubitPoint **close_point,
167 : : int &num_close,
168 : : int max_close,
169 : : int min_close );
170 : : //- return an array of points that are close to the given point
171 : :
172 : : static CubitStatus eval_bezier_patch( CubitFacet *facet, CubitVector &areacoord,
173 : : CubitVector &pt );
174 : : static CubitStatus eval_bezier_patch_normal( CubitFacet *facet,
175 : : CubitVector &areacoord,
176 : : CubitVector &normal );
177 : : static CubitStatus hodograph( CubitFacet *facet,
178 : : CubitVector &areacoord,
179 : : CubitVector Nijk[10] );
180 : : static CubitStatus project_to_patch( CubitFacet *facet,
181 : : CubitVector &ac,
182 : : const CubitVector &pt,
183 : : CubitVector &eval_pt,
184 : : CubitVector *eval_norm,
185 : : CubitBoolean &outside,
186 : : double compare_tol,
187 : : int edge_id /* = -1*/);
188 : : //- Project a point to a bezier patch. Pass in the area
189 : : //- of the point projected to the linear facet. Function
190 : : //- assumes that the point is contained within the patch -
191 : : //- if not, it will project to one of its edges.
192 : :
193 : : static void ac_at_edge( CubitVector &fac,
194 : : CubitVector &eac,
195 : : int edge_id );
196 : : //- determine the area coordinate of the facet at the edge
197 : :
198 : : static CubitBoolean is_at_vertex( CubitFacet *facet,
199 : : const CubitVector &pt,
200 : : CubitVector &ac,
201 : : double compare_tol,
202 : : CubitVector &eval_pt,
203 : : CubitVector *eval_norm_ptr );
204 : :
205 : : static CubitBoolean move_ac_inside( CubitVector &ac, double tol );
206 : : //- find the closest area coordinate to the boundary of the
207 : : //- patch if any of its components are < 0
208 : : //- Return if the ac was modified.
209 : :
210 : : void check_faceting();
211 : :
212 : :
213 : : // For a given point, find which pairs of edges
214 : : // should be C1 across that point.
215 : : CubitStatus mark_edge_pairs(CubitPoint* point);
216 : :
217 : : //loop over the points and call the above function
218 : : CubitStatus pair_edges();
219 : :
220 : : public:
221 : :
222 : : FacetEvalTool();
223 : : //constructor
224 : : ~FacetEvalTool();
225 : : //destructor
226 : :
227 : : CubitStatus initialize( DLIList<CubitFacet*> &facet_list, DLIList<CubitPoint*> &point_list,
228 : : int interp_order/* = -1*/, double min_dot/* = 0.707106781185*/ );
229 : :
230 : 0 : int tool_id()
231 : 0 : { return toolID; };
232 : :
233 : : CubitStatus reverse_facets();
234 : : static CubitStatus reverse_facets(DLIList<CubitFacet *> &facets);
235 : :
236 : 1474 : int get_output_id() { return output_id; }
237 : 275 : void set_output_id( int id ) { output_id = id; }
238 : :
239 : : CubitStatus save( FILE *fp );
240 : : // save to a cubit file
241 : : CubitStatus restore( FILE *fp, unsigned int endian,
242 : : int num_facets, int num_edges,
243 : : int num_points, CubitFacet **facets,
244 : : CubitFacetEdge **edges, CubitPoint **points );
245 : :
246 : : //- get and set the interpolation order
247 : : CubitBox bounding_box();
248 : : void set_bounding_box( CubitBox &box ) { *myBBox = box; }
249 : : void reset_bounding_box();
250 : : //- Returns the bounding box for the set of facets (based on the points
251 : : //- used in the faceted set.
252 : :
253 : : double area();
254 : : //- return the surface area of the facets
255 : :
256 : : void calculate_area();
257 : : //- (re) calcualte area to make sure it is up to date.
258 : :
259 : 0 : double get_min_dot() {return minDot;};
260 : : //- return the mindot (feature angle)
261 : :
262 : : CubitStatus closest_point(CubitVector &this_point,
263 : : CubitVector *closest_point_ptr,
264 : : CubitVector *normal_ptr = NULL);
265 : : //- Finds the closest point from the vector (this_point) to the
266 : : //- set of facets that lies on the set of facets. If the point
267 : : //- lies outside this set, the closest point will be on the plane
268 : : //- of the closest facet. The closest_point is set to be that point.
269 : :
270 : : CubitStatus closest_point_trimmed( CubitVector &this_point,
271 : : CubitVector *closest_point_ptr,
272 : : CubitBoolean &lies_outside,
273 : : CubitVector *normal_ptr = NULL);
274 : : //- Finds the closest point from the vector (this_point) to the
275 : : //- set of facets that lies on the set of facets. If the point
276 : : //- lies outside this set, the closest point will be on the edge
277 : : //- of the closest facet.
278 : : //- This function also determines if the point is outside or
279 : : //- inside the current set of facets. You should call
280 : : //- a bounding box test first before this...
281 : :
282 : : CubitFacet* closest_facet( const CubitVector& point );
283 : :
284 : : int is_flat();
285 : : //- Determine if the set of facets are flat (all in the same plane)
286 : :
287 : : static CubitStatus eval_facet( CubitFacet *facet,
288 : : CubitVector &areacoord,
289 : : CubitVector *eval_point,
290 : : CubitVector *eval_normal );
291 : : //- evaluate the facet at the areacoords
292 : :
293 : : CubitStatus eval_facet( CubitFacet *facet,
294 : : CubitVector &pt,
295 : : CubitVector &areacoord,
296 : : CubitVector &close_point,
297 : : CubitBoolean &outside_facet );
298 : : //- evaluate a location and normal on (or near) a facet based on the
299 : : //- the facet area coordinates
300 : :
301 : : static CubitStatus project_to_facet( CubitFacet *facet,
302 : : const CubitVector &pt,
303 : : CubitVector &areacoord,
304 : : CubitVector &close_point,
305 : : CubitBoolean &outside,
306 : : double compare_tol);
307 : : //- project a point to the given facet (assumes that it contains
308 : : //- the point) areacoord is an initial guess - may be changed.
309 : :
310 : : static CubitStatus project_to_facets(DLIList <CubitFacet *> &facet_list,
311 : : CubitFacet *&last_facet,
312 : : int interp_order,
313 : : double compare_tol,
314 : : const CubitVector &this_point,
315 : : CubitBoolean trim,
316 : : CubitBoolean *outside,
317 : : CubitVector *closest_point_ptr,
318 : : CubitVector *normal_ptr);
319 : : //- project a point to the a list of facets using the interpolation
320 : : //- order. Option to trim to boundary
321 : :
322 : : static CubitStatus eval_edge( CubitFacet *facet,
323 : : int vert0, int vert1,
324 : : CubitVector &pt_on_plane,
325 : : CubitVector &close_point);
326 : : static CubitStatus project_to_facetedge( CubitFacet *facet,
327 : : int vert0, int vert1,
328 : : const CubitVector &the_point,
329 : : CubitVector &pt_on_plane,
330 : : CubitVector &close_point,
331 : : CubitBoolean &outside_facet,
332 : : CubitBoolean must_be_on_edge = CUBIT_FALSE);
333 : : //- project a point to the edge of a facet
334 : : //- if must_be_on_edge is true, close_point will be on the facet edge
335 : :
336 : : static CubitStatus eval_point( CubitFacet *facet,
337 : : int vertex_id,
338 : : CubitVector &close_point );
339 : : //- evaluate a location at a vertex of a facet
340 : :
341 : : static CubitStatus eval_facet_normal( CubitFacet *facet,
342 : : CubitVector &areacoord,
343 : : CubitVector &normal );
344 : : //- evaluate the normal on a facet (use the interpOrder)
345 : :
346 : : static void project_to_facet_plane( CubitFacet *facet,
347 : : const CubitVector &pt,
348 : : CubitVector &point_on_plane,
349 : : double &dist );
350 : : //- project a point to the plane of a facet
351 : :
352 : : static void facet_area_coordinate( CubitFacet *facet,
353 : : const CubitVector &pt_on_plane,
354 : : CubitVector &areacoord );
355 : : //- define the area coordinates of a point on a plane of the facet
356 : :
357 : : static CubitBoolean is_outside( CubitFacet *facet,
358 : : CubitVector &areacoord );
359 : : //- determines if the point at the areacoord is outside
360 : : //- the range of facets
361 : :
362 : 0 : DLIList<DLIList<CubitFacetEdge*>*> *loops() { return &myLoopList; };
363 : : //- return the edge loop list
364 : :
365 : 1584 : int interp_order() { return interpOrder; };
366 : : //- return the interpolation order of the facet eval tool
367 : :
368 : 1232 : void get_facets( DLIList<CubitFacet*> &facet_list ) { facet_list += myFacetList; };
369 : 0 : void get_edges( DLIList<CubitFacetEdge*> &edge_list) { edge_list += myEdgeList; };
370 : 1364 : void get_points( DLIList<CubitPoint*> &point_list ) { point_list += myPointList; };
371 : : void remove_facets( DLIList<CubitFacet*>& facet_list );
372 : : //- return the facets and points defining the surface (append to list)
373 : :
374 : : void replace_facets( DLIList<CubitFacet *> &facet_list );
375 : : // replace the facets on this FacetEvalTool
376 : :
377 : : int num_facets() const {return myFacetList.size();};
378 : : //- get the number of facets
379 : :
380 : : int num_points() const {return myPointList.size();};
381 : : //- get the number of points
382 : :
383 : : double compare_tol() { return compareTol; };
384 : : void compare_tol(double tol) { compareTol = tol; };
385 : : //- set/get the comparison tolerance.
386 : :
387 : : static double get_grid_search_time()
388 : : { return timeGridSearch;}
389 : : static double get_facet_project_time()
390 : : { return timeFacetProject;}
391 : : //- Returns static timing information for the two functions.
392 : :
393 : : void debug_draw_facets(int color = -1 );
394 : :
395 : : void transform_control_points( CubitTransformMatrix &rotmat );
396 : : // transform all control points on the surface a specified
397 : : // transformation matrix
398 : :
399 : : CubitStatus init_bezier_surface( );
400 : : static CubitStatus init_bezier_facet( CubitFacet *facet );
401 : : static CubitStatus init_bezier_edge( CubitFacetEdge *edge, double min_dot );
402 : : static CubitStatus init_boundary_bezier_edge( CubitFacetEdge *edge );
403 : : static CubitStatus init_edge_control_points( CubitVector &P0,
404 : : CubitVector &P3,
405 : : CubitVector &N0,
406 : : CubitVector &N3,
407 : : CubitVector &T0,
408 : : CubitVector &T3,
409 : : CubitVector Pi[3] );
410 : : static CubitStatus init_edge_control_points_single(
411 : : CubitVector &P0,
412 : : CubitVector &P3,
413 : : CubitVector &T0,
414 : : CubitVector &T3,
415 : : CubitVector Pi[3] );
416 : : static CubitStatus init_facet_control_points( CubitVector N[6],
417 : : CubitVector P[3][5],
418 : : CubitVector G[6] );
419 : : //- Initialize the bezier control points for
420 : : //- G1 bezier patch surface interpolation
421 : :
422 : : CubitStatus parameterize();
423 : : //- define a parameterization
424 : :
425 : : static CubitStatus compute_curve_tangent( CubitFacetEdge *edge,
426 : : double min_dot,
427 : : CubitVector &T0,
428 : : CubitVector &T3 );
429 : : //- compute curve tangent vectors at an aedge on the boundary
430 : : static CubitFacetEdge *next_boundary_edge( CubitFacetEdge *this_edge, CubitPoint *p0 );
431 : : //- return the next edge on the boundary
432 : :
433 : : CubitStatus get_intersections(CubitVector point1,
434 : : CubitVector point2,
435 : : DLIList<CubitVector*>& intersection_list,
436 : : bool bounded = CUBIT_FALSE);
437 : :
438 : :
439 : : // The following copyright applies to the following two functions...
440 : : //
441 : : // Copyright 2001, softSurfer (www.softsurfer.com)
442 : : //
443 : : // This code may be freely used and modified for any purpose
444 : : // providing that this copyright notice is included with it.
445 : : // SoftSurfer makes no warranty for this code, and cannot be held
446 : : // liable for any real or imagined damage resulting from its use.
447 : : // Users of this code must verify correctness for their application.
448 : :
449 : : static int intersect_ray( CubitVector &origin, CubitVector &direction, CubitFacet* facet, CubitVector* point, double &distance );
450 : : //- Find intersection point of a ray and a facet
451 : : // Return: -1 = triangle is degenerate (a segment or point)
452 : : // 0 = disjoint (no intersect)
453 : : // 1 = intersect at unique point
454 : : // 2 = are in the same plane
455 : :
456 : : static int intersect_ray( CubitVector &origin, CubitVector &direction, CubitFacetEdge* facet, CubitVector* point, double &distance );
457 : : //- Find intersection point of a ray and a facet edge
458 : : // Return: -1 = edge is degenerate (a point)
459 : : // 0 = disjoint (no intersect)
460 : : // 1 = intersect at unique point
461 : : // 2 = are the same line (infinite points)
462 : :
463 : :
464 : : CubitStatus get_loops_from_facets(DLIList<CubitFacetEdge*> &all_edge_list,
465 : : DLIList<DLIList<CubitFacetEdge*>*> &loop_list );
466 : :
467 : :
468 : : static double contained_volume( DLIList<CubitFacet *> &facets );
469 : : //- computed the contained volume of a set of facets
470 : : // Return: volume of facets (may be negative based on orientation of facets)
471 : :
472 : : };
473 : :
474 : : #endif // SMOOTH_FACET_EVAL_TOOL_HPP
475 : :
476 : :
477 : :
478 : :
479 : :
|
