cgma
FacetEvalTool.hpp
Go to the documentation of this file.
00001 //- Class: FacetEvalTool
00002 //- Description:  The FacetEvalTool is a general purpose tool that uses a set
00003 //-               of facets to determine certain geometry computations.  This
00004 //-               class can be used to define a surface or do fast bounding box
00005 //-               calculations, or determine point inside or outside.
00006 //- Assumptions:  It is assumed/required that the facets be continuous through
00007 //-               the set.  In other words they must shared nodes and be conected
00008 //-               through each other.  The algorithms assume, for efficient searching,
00009 //-               that from one facet, you can get to any other facet in the set.
00010 //- Owner: Steven J. Owen
00011 //- Checked by: 
00012 
00013 #ifndef SMOOTH_FACET_EVAL_TOOL_HPP
00014 #define SMOOTH_FACET_EVAL_TOOL_HPP
00015 
00016 #include "DLIList.hpp"
00017 #include "AbstractTree.hpp"
00018 #include "CubitFacet.hpp" //For inline function.
00019 
00020 template <class Y> class KDDTree;
00021 
00022 namespace FacetEvalToolUtils {
00023 
00024 template <typename R> inline R determ3(R p1, R q1, R p2, R q2, R p3, R q3) {
00025   return (q3) * ((p2) - (p1)) + (q2) * ((p1) - (p3)) + (q1) * ((p3) - (p2));
00026 }
00027 template <typename R> inline R sqr(R a) { return a * a; }
00028 template <typename R> inline R cube(R a) { return sqr(a) * a; }
00029 template <typename R> inline R quart(R a) { return sqr(a) * sqr(a); }
00030 template <typename R> inline R blend(R x) {
00031   return -2.0 * (x) * (x) * (x) + 3.0 * (x) * (x);
00032 }
00033 
00034 } // namespace FacetEvalToolUtils
00035 
00036 class CubitBox;
00037 class CubitFacet;
00038 class CubitPoint;
00039 class CubitVector;
00040 class CubitTransformMatrix;
00041 
00042 class FacetEvalTool
00043 {
00044 private:
00045   static double timeGridSearch;
00046   static double timeFacetProject;
00047   static int numEvals;
00048 
00049   int toolID;
00050   int interpOrder;
00051     //- interpolation order 0=linear, 1=gradient, 2=quadratic, 4=quartic Bezier
00052   DLIList<CubitFacet*> myFacetList;
00053     //- Facets that are used for evaluation.
00054   DLIList<CubitFacet*> markedFacetList;
00055     //- Facets that are used for evaluation.
00056   DLIList<CubitPoint*> myPointList;
00057     //- point list
00058   DLIList<CubitFacetEdge*> myEdgeList;
00059   DLIList<DLIList<CubitFacetEdge*>*> myLoopList;
00060     //- edges only used for bezier interp
00061   AbstractTree <CubitFacet*> *aTree;
00062     //- used for surfaces with many facets
00063   CubitBox *myBBox;
00064     //- Bounding box of facets
00065   double myArea;
00066     //- area of the facet surface
00067   double compareTol;
00068     //- for comparing bounding box
00069   CubitFacet *lastFacet;
00070     //- The last facet evaluated
00071   int isFlat;
00072     //- The surface represented by the facets is flat
00073   CubitBoolean isParameterized;
00074     //- Surface contains a parameterization (uv values initialized)
00075   double minDot;  
00076     //- computed from feature angle
00077   int output_id;
00078 
00079   bool have_data_to_calculate_bbox(void);
00080 
00081   void set_up_grid_search(double geom_factor);
00089   void facets_from_search_grid( CubitVector &this_point,
00090                                 DLIList<CubitFacet *> &facet_list,
00091                                 double &tol_used );
00099   void facets_from_search_grid( CubitVector &this_point,
00100                                 double compare_tol,
00101                                 DLIList<CubitFacet *> &facet_list );
00102     //- set up a grid search data structure if we have lots of facets
00103 
00104   CubitStatus get_points_from_facets(DLIList<CubitFacet*> &facet_list,
00105                                      DLIList<CubitPoint*> &point_list );
00106     //- populates the point_list with points contained by the given facets.
00107   CubitStatus get_edges_from_facets(DLIList<CubitFacet*> &facet_list,
00108                                     DLIList<CubitFacetEdge*> &edge_list );
00109     //- populates the edge_list with edges contained by the given facets
00110 
00111     //- populates the loop_list with edges contained by the given edge list
00112   void destroy_facets();
00113     //- Destroys the facets, and points.
00114 
00115   void debug_draw_edges(int color = -1 );
00116   void debug_draw_edge(CubitFacetEdge *edge, int color = -1 );
00117   void debug_draw_vec( CubitVector & vec, int color = -1 );
00118   void debug_draw_facet_normals( int color = -1 );
00119   void debug_draw_point_normals( int color = -1 );
00120   void debug_draw_bezier_edges( int color = -1 );
00121   void debug_draw_bezier_facets( int color = -1 );
00122   void debug_draw_bezier_facet( CubitFacet *facet, int color = -1 );
00123   void debug_draw_line( CubitVector &begin, CubitVector &end, int color = -1 );
00124   void debug_draw_eval_bezier_facet( CubitFacet *facet );
00125   void debug_draw_location( CubitVector &location, int color = -1 );
00126   void write_loops();
00127     //- drawing functions for debugging.
00128 
00129   CubitStatus project_to_facets(CubitVector &this_point,
00130                                 CubitBoolean trim,
00131                                 CubitBoolean *outside,
00132                                 CubitVector *closest_point_ptr,
00133                                 CubitVector *normal_ptr);
00134 
00135   CubitStatus init_gradient();
00136     //- initiaize the gradients (tangent planes at the points)
00137 
00138   CubitStatus init_quadrics();
00139     //- initialize quadric coefficients at the points
00140 
00141   CubitStatus eval_quadratic( CubitFacet *facet, 
00142                               int pt_idx, 
00143                               CubitVector &eval_pt,
00144                               CubitVector &qpoint,
00145                               CubitVector &qnorm );
00146     //- evaluate on facet based on quadratic approximation
00147 
00148   void on_circle( CubitVector &ptA,
00149                   CubitVector &normA,
00150                   CubitVector &ptB,
00151                   CubitVector &normB,
00152                   CubitVector &eval_pt,
00153                   CubitVector &pt_on_circle,
00154                   CubitVector &normal );
00155     //- compute a projection of a point onto a circle.  The circle
00156     //- is defined by two points and two normals.  Return the 
00157     //- projected point and its normal 
00158 
00159   void eval_quadric( CubitFacet *facet,
00160                      int pt_index,
00161                      CubitVector &eval_pt,
00162                      CubitVector &qpt );
00163     //- evaluate on facet based on quadratic approximation
00164 
00165   CubitStatus get_close_points( CubitPoint *point,
00166                                 CubitPoint **close_point,
00167                                 int &num_close,
00168                                 int max_close,
00169                                 int min_close );
00170     //- return an array of points that are close to the given point
00171 
00172   static CubitStatus eval_bezier_patch( CubitFacet *facet, CubitVector &areacoord,
00173                                CubitVector &pt );
00174   static CubitStatus eval_bezier_patch_normal( CubitFacet *facet, 
00175                                                CubitVector &areacoord,
00176                                                CubitVector &normal );
00177   static CubitStatus hodograph( CubitFacet *facet, 
00178                                 CubitVector &areacoord,
00179                                 CubitVector Nijk[10] );
00180   static CubitStatus project_to_patch( CubitFacet *facet,                                                
00181                                 CubitVector &ac,
00182                                 const CubitVector &pt,
00183                                 CubitVector &eval_pt,
00184                                 CubitVector *eval_norm,
00185                                 CubitBoolean &outside,
00186                                 double compare_tol,
00187                                 int edge_id /* = -1*/);
00188     //- Project a point to a bezier patch. Pass in the area
00189     //- of the point projected to the linear facet.  Function
00190     //- assumes that the point is contained within the patch -
00191     //- if not, it will project to one of its edges.
00192 
00193   static void ac_at_edge( CubitVector &fac,
00194                           CubitVector &eac,
00195                           int edge_id );
00196     //- determine the area coordinate of the facet at the edge
00197 
00198   static CubitBoolean is_at_vertex( CubitFacet *facet,
00199                                     const CubitVector &pt,
00200                                     CubitVector &ac,
00201                                     double compare_tol,
00202                                     CubitVector &eval_pt,
00203                                     CubitVector *eval_norm_ptr );
00204 
00205   static CubitBoolean move_ac_inside( CubitVector &ac, double tol );
00206     //-  find the closest area coordinate to the boundary of the
00207     //-  patch if any of its components are < 0
00208     //-  Return if the ac was modified.
00209 
00210   void check_faceting();
00211 
00212   
00213     // For a given point, find which pairs of edges
00214     // should be C1 across that point.
00215   CubitStatus mark_edge_pairs(CubitPoint* point);
00216 
00217     //loop over the points and call the above function
00218   CubitStatus pair_edges();
00219   
00220 public:
00221   
00222   FacetEvalTool();
00223   //constructor
00224   ~FacetEvalTool();
00225     //destructor
00226 
00227   CubitStatus initialize( DLIList<CubitFacet*> &facet_list, DLIList<CubitPoint*> &point_list,
00228                           int interp_order/* = -1*/, double min_dot/* = 0.707106781185*/ );
00229 
00230   int tool_id()
00231     { return toolID; };
00232 
00233   CubitStatus reverse_facets();
00234   static CubitStatus reverse_facets(DLIList<CubitFacet *> &facets);
00235   
00236   int get_output_id() { return output_id; }
00237   void set_output_id( int id ) { output_id = id; }
00238 
00239   CubitStatus save( FILE *fp );
00240     // save to a cubit file
00241   CubitStatus restore( FILE *fp, unsigned int endian,
00242                        int num_facets, int num_edges, 
00243                        int num_points, CubitFacet **facets, 
00244                        CubitFacetEdge **edges, CubitPoint **points );
00245   
00246     //- get and set the interpolation order
00247   CubitBox bounding_box();
00248   void set_bounding_box( CubitBox &box ) { *myBBox = box; }
00249   void reset_bounding_box();
00250     //- Returns the bounding box for the set of facets (based on the points
00251     //- used in the faceted set.
00252 
00253   double area();
00254     //- return the surface area of the facets
00255 
00256   void calculate_area();
00257     //- (re) calcualte area to make sure it is up to date.
00258 
00259   double get_min_dot() {return minDot;};
00260     //- return the mindot (feature angle)
00261 
00262   CubitStatus closest_point(CubitVector &this_point, 
00263                             CubitVector *closest_point_ptr,
00264                             CubitVector *normal_ptr = NULL);
00265     //- Finds the closest point from the vector (this_point) to the
00266     //- set of facets that lies on the set of facets.  If the point
00267     //- lies outside this set, the closest point will be on the plane
00268     //- of the closest facet.  The closest_point is set to be that point.
00269   
00270   CubitStatus closest_point_trimmed( CubitVector &this_point, 
00271                                      CubitVector *closest_point_ptr,
00272                                      CubitBoolean &lies_outside,
00273                                      CubitVector *normal_ptr = NULL);
00274     //- Finds the closest point from the vector (this_point) to the
00275     //- set of facets that lies on the set of facets.  If the point
00276     //- lies outside this set, the closest point will be on the edge
00277     //- of the closest facet.  
00278     //- This function also determines if the point is outside or
00279     //- inside the current set of facets.  You should call
00280     //- a bounding box test first before this...
00281 
00282   CubitFacet* closest_facet( const CubitVector& point );
00283 
00284   int is_flat();
00285     //- Determine if the set of facets are flat (all in the same plane)
00286 
00287   static CubitStatus eval_facet( CubitFacet *facet, 
00288                                  CubitVector &areacoord,
00289                                  CubitVector *eval_point,
00290                                  CubitVector *eval_normal );
00291     //- evaluate the facet at the areacoords
00292 
00293   CubitStatus eval_facet( CubitFacet *facet, 
00294                           CubitVector &pt,
00295                           CubitVector &areacoord, 
00296                           CubitVector &close_point,
00297                           CubitBoolean &outside_facet );
00298     //- evaluate a location and normal on (or near) a facet based on the
00299     //- the facet area coordinates
00300 
00301   static CubitStatus project_to_facet( CubitFacet *facet, 
00302                                        const CubitVector &pt,
00303                                        CubitVector &areacoord,
00304                                        CubitVector &close_point,
00305                                        CubitBoolean &outside,
00306                                        double compare_tol);
00307     //- project a point to the given facet (assumes that it contains
00308     //- the point) areacoord is an initial guess - may be changed.
00309 
00310   static CubitStatus project_to_facets(DLIList <CubitFacet *> &facet_list,
00311                                  CubitFacet *&last_facet,
00312                                  int interp_order,
00313                                  double compare_tol,
00314                                  const CubitVector &this_point,
00315                                  CubitBoolean trim,
00316                                  CubitBoolean *outside,
00317                                  CubitVector *closest_point_ptr,
00318                                  CubitVector *normal_ptr);
00319     //- project a point to the a list of facets using the interpolation
00320     //- order.  Option to trim to boundary
00321 
00322   static CubitStatus eval_edge( CubitFacet *facet, 
00323                          int vert0, int vert1,
00324                          CubitVector &pt_on_plane, 
00325                          CubitVector &close_point);
00326   static CubitStatus project_to_facetedge( CubitFacet *facet, 
00327                          int vert0, int vert1,
00328                          const CubitVector &the_point,
00329                          CubitVector &pt_on_plane, 
00330                          CubitVector &close_point,
00331                          CubitBoolean &outside_facet, 
00332                          CubitBoolean must_be_on_edge = CUBIT_FALSE);
00333     //- project a point to the edge of a facet
00334     //- if must_be_on_edge is true, close_point will be on the facet edge
00335 
00336   static CubitStatus eval_point( CubitFacet *facet, 
00337                           int vertex_id, 
00338                           CubitVector &close_point );
00339     //- evaluate a location at a vertex of a facet
00340 
00341   static CubitStatus eval_facet_normal( CubitFacet *facet,
00342                                  CubitVector &areacoord,
00343                                  CubitVector &normal );
00344     //- evaluate the normal on a facet (use the interpOrder)
00345 
00346   static void project_to_facet_plane( CubitFacet *facet,
00347                                const CubitVector &pt,
00348                                CubitVector &point_on_plane,
00349                                double &dist );
00350     //- project a point to the plane of a facet
00351 
00352   static void facet_area_coordinate( CubitFacet *facet,
00353                                      const CubitVector &pt_on_plane,
00354                                      CubitVector &areacoord );
00355     //- define the area coordinates of a point on a plane of the facet
00356 
00357   static CubitBoolean is_outside( CubitFacet *facet, 
00358                            CubitVector &areacoord );
00359     //- determines if the point at the areacoord is outside 
00360     //- the range of facets
00361 
00362   DLIList<DLIList<CubitFacetEdge*>*> *loops() { return &myLoopList; };
00363     //- return the edge loop list
00364 
00365   int interp_order() { return interpOrder; };
00366     //- return the interpolation order of the facet eval tool
00367 
00368   void get_facets( DLIList<CubitFacet*> &facet_list ) { facet_list += myFacetList; };
00369   void get_edges( DLIList<CubitFacetEdge*> &edge_list) { edge_list += myEdgeList; };
00370   void get_points( DLIList<CubitPoint*> &point_list ) { point_list += myPointList; };
00371   void remove_facets( DLIList<CubitFacet*>& facet_list );
00372     //- return the facets and points defining the surface (append to list)
00373 
00374   void replace_facets( DLIList<CubitFacet *> &facet_list );
00375     // replace the facets on this FacetEvalTool
00376 
00377   int num_facets() const {return myFacetList.size();};
00378     //- get the number of facets
00379       
00380   int num_points() const {return myPointList.size();};
00381     //- get the number of points
00382 
00383   double compare_tol() { return compareTol; };
00384   void compare_tol(double tol) { compareTol = tol; };
00385     //- set/get the comparison tolerance.
00386 
00387   static double get_grid_search_time()
00388     { return timeGridSearch;}
00389   static double get_facet_project_time()
00390     { return timeFacetProject;}
00391     //- Returns static timing information for the two functions.
00392 
00393   void debug_draw_facets(int color = -1 );
00394 
00395   void transform_control_points( CubitTransformMatrix &rotmat );
00396     // transform all control points on the surface a specified 
00397     // transformation matrix
00398 
00399   CubitStatus init_bezier_surface( );
00400   static CubitStatus init_bezier_facet( CubitFacet *facet );
00401   static CubitStatus init_bezier_edge( CubitFacetEdge *edge, double min_dot );
00402   static CubitStatus init_boundary_bezier_edge( CubitFacetEdge *edge );
00403   static CubitStatus init_edge_control_points( CubitVector &P0,
00404                                                CubitVector &P3,
00405                                                CubitVector &N0,
00406                                                CubitVector &N3,
00407                                                CubitVector &T0,
00408                                                CubitVector &T3,
00409                                                CubitVector Pi[3] );
00410   static CubitStatus init_edge_control_points_single( 
00411                               CubitVector &P0,
00412                               CubitVector &P3,
00413                               CubitVector &T0,
00414                               CubitVector &T3,
00415                               CubitVector Pi[3] );                                      
00416   static CubitStatus init_facet_control_points( CubitVector N[6],
00417                                                 CubitVector P[3][5],
00418                                                 CubitVector G[6] );
00419     //- Initialize the bezier control points for 
00420     //- G1 bezier patch surface interpolation
00421 
00422   CubitStatus parameterize();
00423     //- define a parameterization
00424 
00425   static CubitStatus compute_curve_tangent( CubitFacetEdge *edge,
00426                                      double min_dot,
00427                                      CubitVector &T0,
00428                                      CubitVector &T3 );
00429     //- compute curve tangent vectors at an aedge on the boundary
00430   static CubitFacetEdge *next_boundary_edge( CubitFacetEdge *this_edge, CubitPoint *p0 );
00431     //- return the next edge on the boundary
00432 
00433   CubitStatus get_intersections(CubitVector point1,
00434                                 CubitVector point2,
00435                                 DLIList<CubitVector*>& intersection_list,
00436                                 bool bounded = CUBIT_FALSE);
00437 
00438 
00439   // The following copyright applies to the following two functions...
00440   //
00441   // Copyright 2001, softSurfer (www.softsurfer.com)
00442   //
00443   // This code may be freely used and modified for any purpose
00444   // providing that this copyright notice is included with it.
00445   // SoftSurfer makes no warranty for this code, and cannot be held
00446   // liable for any real or imagined damage resulting from its use.
00447   // Users of this code must verify correctness for their application.
00448 
00449   static int intersect_ray( CubitVector &origin, CubitVector &direction, CubitFacet* facet, CubitVector* point, double &distance );
00450     //- Find intersection point of a ray and a facet
00451     //    Return: -1 = triangle is degenerate (a segment or point)
00452     //             0 = disjoint (no intersect)
00453     //             1 = intersect at unique point
00454     //             2 = are in the same plane
00455 
00456   static int intersect_ray( CubitVector &origin, CubitVector &direction, CubitFacetEdge* facet, CubitVector* point, double &distance );
00457     //- Find intersection point of a ray and a facet edge
00458     //    Return: -1 = edge is degenerate (a point)
00459     //             0 = disjoint (no intersect)
00460     //             1 = intersect at unique point
00461     //             2 = are the same line (infinite points)
00462 
00463   
00464   CubitStatus get_loops_from_facets(DLIList<CubitFacetEdge*> &all_edge_list,
00465                                     DLIList<DLIList<CubitFacetEdge*>*> &loop_list );
00466   
00467   
00468   static double contained_volume( DLIList<CubitFacet *> &facets );
00469     //- computed the contained volume of a set of facets
00470     //    Return: volume of facets (may be negative based on orientation of facets)
00471   
00472 };
00473 
00474 #endif // SMOOTH_FACET_EVAL_TOOL_HPP
00475 
00476 
00477 
00478     
00479 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines