cgma
AnalyticGeometryTool.hpp
Go to the documentation of this file.
```00001 //-------------------------------------------------------------------------
00003 //
00004 // Filename      : AnalyticGeometryTool.hpp
00005 //
00006 // Purpose       : This class performs calculations on analytic geometry
00007 //                 (points, lines, arcs, planes, polygons).  Capabilities
00008 //                 include vector and point math, matrix operations,
00009 //                 measurements, intersections and comparison/containment checks.
00010 //
00011 // Special Notes : As most of these functions were taken from Caterpillar's
00012 //                 Cmath library with minimal modification for this class,
00013 //                 the vectors, planes, etc. are represented in non-CUBIT
00014 //                 format.  These could be converted with some effort.
00015 //
00016 //                 Other Notes:
00017 //                 ------------
00018 //                 1. Points and vectors are represented as double[3].
00019 //
00020 //                 2. Lines (unless denoted as 'segments') are defined as
00021 //                    a point and a vector (point-direction form).  The parametric
00022 //                    equations of a line are:
00023 //
00024 //                         x = xo + t*a
00025 //                         y = yo + t*a
00026 //                         z = zo + t*a
00027 //
00028 //                 3. Planes are also defined as a point and a vector
00029 //                    (point-normal form).  The equation of a plane is:
00030 //
00031 //                         a*(x-xo) + b*(y-yo) + c*(z-zo) = 0
00032 //
00033 //                     For example, the equation of a plane passing through the
00034 //                     point (3,-1,7) and perpendicular to the vector n = (4,2,-5)
00035 //                     is: 4(x-3) + 2(y+1) -5(z-7) = 0
00036 //
00037 //                     This can be reduced to:
00038 //
00039 //                         A*x + B*y + C*z + D = 0
00040 //
00041 //                     For the example, the equation is: 4x + 2y -5z + 25 = 0
00042 //
00043 //                 4. Most functions which can accept a 3x3 matrix have
00044 //                    an analogous 4x4 function.  This is for convenience
00045 //                    only (the extra rows/columns are ignored within the function).
00046 //
00047 //
00048 // Creator       : Steve Storm
00049 //
00050 // Creation Date : 10/16/98
00051 //-------------------------------------------------------------------------
00052
00053 #ifndef ANALYTIC_GEOMETRY_TOOL_HPP
00054 #define ANALYTIC_GEOMETRY_TOOL_HPP
00055
00056 #include <math.h>
00057 #include "CubitDefines.h"
00058 #include "CGMGeomConfigure.h"
00059
00060 #include <iostream>
00061
00062 class CubitBox;
00063 class CubitPlane;
00064 class CubitVector;
00065 template <class X> class DLIList;
00066
00067   // Multiplication constants
00068   const double AGT_E =             2.7182818284590452354;   // e
00069   const double AGT_LOG_2E =        1.4426950408889634074;   // log2(e)
00070   const double AGT_LOG_10E =       0.43429448190325182765; // log10(e)
00071   const double AGT_LN_2 =          0.69314718055994530942; // ln(2)
00072   const double AGT_LN_10 =         2.30258509299404568402; // ln(10)
00073   const double AGT_PI  =           3.14159265358979323846; // PI
00074   const double AGT_PI_DIV_2 =      1.57079632679489661923; // PI/2
00075   const double AGT_PI_DIV_4 =      0.78539816339744830962; // PI/4
00076   const double AGT_1_DIV_PI =      0.31830988618379067154; // 1/PI
00077   const double AGT_2_DIV_PI =      0.63661977236758134308; // 2/PI
00078   const double AGT_2_DIV_SQRT_PI = 1.12837916709551257390; // 2/(sqrt(PI))
00079   const double AGT_SQRT_2 =        1.41421356237309504880; // sqrt(2)
00080   const double AGT_SQRT_1_DIV_2 =  0.70710678118654752440; // sqrt(1/2)
00081   const double AGT_PI_DIV_180 =    0.017453292519943295;   // PI/180
00082   const double AGT_180_DIV_PI =    57.295779513082323;     // 180/PI
00083
00084   const double DEGtoRAD = AGT_PI_DIV_180;
00085   const double RADtoDEG = AGT_180_DIV_PI;
00086
00087   typedef struct {
00088      double Vec1[3];    /* First vector that defines plane of arc */
00089      double Vec2[3];    /* Second vector that defines plane of arc */
00090      double Origin[3];  /* Center of arc (on plane) */
00091      double StartAngle; /* Starting angle (in radians) of arc */
00092      double EndAngle;   /* End angle (in radians) of arc */
00094   } AGT_3D_Arc;
00095
00096   // Return values */
00097   enum AgtEquality {
00098      AGT_EQUAL       = 0, // define for comparisons
00099      AGT_LESSTHAN    = 1, // define for comparisons
00100      AGT_GREATERTHAN = 2  // define for comparisons
00101   };
00102
00103   enum AgtSide {
00104      AGT_ON_PLANE =  0,  // define for which side of plane
00105      AGT_POS_SIDE =  1,  // define for which side of plane
00106      AGT_NEG_SIDE =  2,  // define for which side of plane
00107      AGT_INT_PLANE = 3,  // define for intersects plane
00108      AGT_CROSS = 4       // define for line crossing plane
00109   };
00110
00111   enum AgtOrientation {
00112     AGT_OPP_DIR = 0,  // define for vector direction comparision
00113     AGT_SAME_DIR = 1  // define for vector direction comparision
00114   };
00115
00116   enum AgtDistanceMethod {
00117     AGT_FRACTION = 0,  // define for distance along a curve
00118     AGT_ABSOLUTE = 1   // define for distance along a curve
00119   };
00120
00121 // Macros
00122 #define number_sign(number) (number<0 ? -1 : 1)
00123 #define copy_vec(vec1,vec2) copy_pnt(vec1,vec2)
00124 #define set_vec(pnt,x,y,z) set_pnt(pnt,x,y,z)
00125
00127 //                               MAGIC SOFTWARE                              //
00128 // See note later in file regarding MAGIC softare.                           //
00130 // Minimal 2D bounding box - parameters
00131 const int maxPartition = 32;
00132 const int invInterp = 32;
00133 const double angleMin = -0.25*AGT_PI;
00134 const double angleMax = +0.25*AGT_PI;
00135 const double angleRange = 0.5*AGT_PI;
00136
00137 typedef struct
00138 {
00139     double x, y;
00140 }
00141 Point2;
00142
00143 typedef struct
00144 {
00145     double x, y, z;
00146 }
00147 Point3;
00148
00149 typedef struct
00150 {
00151     Point2 center;
00152     Point2 axis[2];
00153     double extent[2];
00154 }
00155 OBBox2;
00156
00157 typedef struct
00158 {
00159     Point3 center;
00160     Point3 axis[3];
00161     double extent[3];
00162 }
00163 OBBox3;
00164
00165 // threshold on determinant to conclude lines/rays/segments are parallel
00166 const double par_tolerance = 1e-06;
00167 #define DIST(x) (Sqrt(Abs(x)))
00168
00169 // Line in 2D
00170 typedef struct
00171 {
00172     // line is Dot(N,X) = c, N not necessarily unit length
00173     Point2 N;
00174     double c;
00175 }
00176 AgtLine;
00177
00178 // Triangle in 2D
00179 typedef struct
00180 {
00181     Point2 v[3];
00182 }
00183 Triangle;
00184
00185 // Linked list of triangle
00186 typedef struct AgtTriList
00187 {
00188     Triangle* tri;
00189     AgtTriList* next;
00190 }
00191 AgtTriList;
00192
00193 //---------------------------------------------------------------------------
00194 // Lines in 3D
00195 //---------------------------------------------------------------------------
00196 typedef struct
00197 {
00198     // Line is L(t) = b+t*m for any real-valued t
00199     // Ray has constraint t >= 0, b is the origin of the ray
00200     // Line segment has constraint 0 <= t <= 1, b and b+m are end points
00201
00202     Point3 b, m;
00203 }
00204 Line3;
00205 //---------------------------------------------------------------------------
00206
00207 //---------------------------------------------------------------------------
00208 // Circles in 3D
00209 //---------------------------------------------------------------------------
00210 typedef struct
00211 {
00212     // Plane containing circle is Dot(N,X-C) = 0 where X is any point in the
00213     // plane.  Vectors U, V, and N form an orthonormal right-handed set
00214     // (matrix [U V N] is orthonormal and has determinant 1).  Circle within
00215     // the plane is parameterized by X = C + R*(cos(A)*U + sin(A)*V) where
00216     // A is an angle in [0,2*pi).
00217
00218     Point3 U, V, N;  // coordinate system of plane containing circle
00219     Point3 C;  // center
00220     double R;  // radius > 0
00221 }
00222 Circle3;
00223 //---------------------------------------------------------------------------
00224
00225 //---------------------------------------------------------------------------
00226 // Triangles in 3D
00227 //---------------------------------------------------------------------------
00228 typedef struct
00229 {
00230     // Triangle points are tri(s,t) = b+s*e0+t*e1 where 0 <= s <= 1,
00231     // 0 <= t <= 1, and 0 <= s+t <= 1.
00232
00233     // If the vertices of the triangle are v0, v1, and v2, then b = v0,
00234     // e0 = v1 - v0, and e1 = v2 - v0.
00235
00236     Point3 b, e0, e1;
00237 }
00238 Triangle3;
00239 //---------------------------------------------------------------------------
00240
00241 //---------------------------------------------------------------------------
00242 // Parallelograms in 3D
00243 //---------------------------------------------------------------------------
00244 typedef struct
00245 {
00246     // Rectoid points are rect(s,t) = b+s*e0+t*e1 where 0 <= s <= 1
00247     // and 0 <= t <= 1.  Could have called this a parallelogram, but
00248     // I do not like typing long words...
00249
00250     Point3 b, e0, e1;
00251 }
00252 Rectoid3;
00253 //---------------------------------------------------------------------------
00254
00255 //---------------------------------------------------------------------------
00256 // Planes in 3D
00257 //---------------------------------------------------------------------------
00258 typedef struct
00259 {
00260     // plane is Dot(N,X) = d
00261     Point3 N;
00262     double d;
00263 }
00264 Plane3;
00265 //---------------------------------------------------------------------------
00266
00267 // END OF MAGIC SOFTWARE
00269
00271 class CUBIT_GEOM_EXPORT AnalyticGeometryTool
00272 {
00273 public:
00274
00275    ~AnalyticGeometryTool();
00276    static AnalyticGeometryTool* instance();
00277
00278    static void delete_instance()
00279    {
00280      if( instance_ )
00281        delete instance_;
00282      instance_ = NULL;
00283    };
00284
00285    //*********************************************************
00286    // Double numbers
00287    //*********************************************************
00288
00292    CubitBoolean dbl_eq( double val_1, double val_2 );
00293
00296    double get_epsilon();
00299    double set_epsilon( double new_epsilon );
00300
00304    void round_near_val( double& val );
00305
00306    //**************************************************************************
00307    // Matrices & Transforms
00308    //**************************************************************************
00309    // Note: For these functions the matrix format is defined as follows:
00310    //
00311    //    Consider the transformation from [x,y,z] to [x',y',z']:
00312    //                             _          _
00313    //                            | x1 y1 z1 0 |
00314    //    [x',y',z',1] = [x,y,z,1]| x2 y2 z2 0 |
00315    //                            | x3 y3 z3 0 |
00316    //                                | ox oy oz 1 |
00317    //                                 -          -
00321    void transform_pnt( double m[4][4], double pin[3], double pout[3] );
00322
00326    void transform_vec( double m3[3][3], double vin[3], double vout[3] );
00327
00331    void transform_vec( double m4[4][4], double vin[3], double vout[3] );
00332
00335    void transform_line( double rot_mtx[4][4], double origin[3], double axis[3] );
00338    void transform_line( double rot_mtx[4][4], CubitVector &origin, CubitVector &axis );
00339
00340
00342
00346    void copy_mtx( double from[3][3],double to[3][3] );
00347    void copy_mtx( double from[4][4], double to[4][4] );
00348    void copy_mtx( double from[4][4], double to[3][3] );
00349    void copy_mtx(double from[3][3], double to[4][4], double *origin = NULL );
00351
00353
00354
00355    void create_rotation_mtx( double theta, double v[3], double mtx3x3[3][3] );
00356    void create_rotation_mtx( double theta, double v[3], double mtx4x4[4][4] );
00358
00362    void add_origin_to_rotation_mtx( double rot_mtx[4][4], double origin[3] );
00363
00365
00366    void identity_mtx( double mtx3x3[3][3] );
00367    void identity_mtx( double mtx4x4[4][4] );
00369
00371
00372
00373
00374
00375
00376
00377
00378
00379    void mtx_to_angs( double mtx3x3[3][3], double &ax, double &ay, double &az );
00380    void mtx_to_angs( double mtx4x4[4][4], double &ax, double &ay, double &az );
00382
00383
00385    // This routine rotates a 3x3 system given a transformation matrix.  The
00386    // matrix can be rotated in place by sending same variable in & out, or a
00387    // new matrix can be created.  In the 4x4 case, the translation portion
00388    // of the matrix is unchanged.
00389    void rotate_system( double mtx[3][3], double sys_in[3][3],
00390                        double sys_out[3][3] );
00391    void rotate_system( double mtx[4][4], double sys_in[4][4],
00392                        double sys_out[4][4] );
00394
00396    double det_mtx( double m[3][3] );
00397
00399    // Multiply matrices together.  If any input is NULL, the identity
00400    // matrix is used in its place.
00401    void mult_mtx( double a[3][3],double b[3][3], double d[3][3] );
00402    void mult_mtx( double a[4][4], double b[4][4], double d[4][4] );
00404
00410    CubitStatus inv_mtx_adj( double mtx[3][3], double inv_mtx[3][3] );
00411
00416    CubitStatus inv_trans_mtx( double transf[4][4], double inv_transf[4][4] );
00417
00419
00420
00421    void vecs_to_mtx( double xvec[3], double yvec[3], double zvec[3],
00422                      double matrix[3][3] );
00423    void vecs_to_mtx( double xvec[3], double yvec[3], double zvec[3],
00424                      double origin[3], double matrix[4][4] );
00426
00428
00429    void print_mtx( double mtx[3][3] );
00430    void print_mtx( double mtx[4][4] );
00432
00433    //**************************************************************************
00434    // 3D Points
00435    //**************************************************************************
00438    void copy_pnt( double pnt_in[3], double pnt_out[3] );
00439
00441
00442    void copy_pnt( double pnt_in[3], CubitVector &cubit_vec );
00443    void copy_pnt( CubitVector &cubit_vec, double pnt_out[3] );
00445
00447    void set_pnt( double pnt[3], double x, double y, double z );
00448
00453    CubitBoolean pnt_eq( double pnt1[3],double pnt2[3] );
00454
00460    CubitStatus mirror_pnt( double pnt[3], double pln_orig[3],
00461                            double pln_norm[3], double m_pnt[3]);
00462
00467    CubitStatus next_pnt( double str_pnt[3], double vec_dir[3], double len,
00468                          double new_pnt[3]);
00469
00470    //***************************************************************************
00471    // 3D Vectors
00472    //***************************************************************************
00474    CubitStatus unit_vec( double vin[3], double vout[3] );
00475
00480    double dot_vec( double uval[3], double vval[3] );
00481
00486    void cross_vec( double uval[3], double vval[3], double cross[3] );
00487
00489    void cross_vec_unit( double uval[3], double vval[3], double cross[3] );
00490
00492    void orth_vecs( double uvect[3], double vvect[3], double wvect[3] );
00493
00496    double mag_vec( double vec[3] );
00497
00500    CubitStatus get_vec ( double str_pnt[3], double stp_pnt[3],
00501                          double vector_out[3] );
00502
00505    CubitStatus get_vec_unit( double str_pnt[3], double stp_pnt[3],
00506                              double uv_out[3] );
00507
00509    void mult_vecxconst( double constant, double vec[3], double vec_out[3] );
00510
00512    void reverse_vec( double vin[3],double vout[3] );
00513
00517    double angle_vec_vec( double v1[3],double v2[3] );
00518
00519    //***************************************************************************
00520    // Distances
00521    //***************************************************************************
00523    double dist_pnt_pnt( double pnt1[3], double pnt2[3] );
00524
00530    double dist_pln_pln( double pln_1_orig[3], double pln_1_norm[3],
00531                         double pln_2_orig[3], double pln_2_norm[3],
00532                         AgtSide *side = NULL, AgtOrientation *orien = NULL,
00533                         unsigned short *status = NULL );
00534
00535    //***************************************************************************
00536    // Intersections
00537    //***************************************************************************
00541    CubitStatus int_ln_pln( double ln_orig[3], double ln_vec[3], double pln_orig[3],
00542                            double pln_norm[3], double int_pnt[3] );
00552    int int_ln_ln( double p1[3], double v1[3], double p2[3], double v2[3],
00553                   double int_pnt1[3], double int_pnt2[3] );
00554
00558    int int_pnt_ln( double pnt[3], double ln_orig[3],  double ln_vec[3],
00559                             double int_pnt[3] );
00560
00564    int int_pnt_pln( double pnt[3], double pln_orig[3], double pln_norm[3],
00565                              double pln_int[3] );
00566
00569    CubitStatus int_pln_pln( double pln_1_orig[3], double pln_1_norm[3],
00570                             double pln_2_orig[3], double pln_2_norm[3],
00571                             double ln_orig[3], double ln_vec[3] );
00572
00573    //**************************************************************************
00574    // Comparison/Containment Tests
00575    //**************************************************************************
00579    CubitBoolean is_vec_par( double vec_1[3], double vec_2[3] );
00580
00583    CubitBoolean is_vec_perp( double vec_1[3],double vec_2[3]);
00584
00587    CubitBoolean is_vecs_same_dir( double vec_1[3], double vec_2[3] );
00588
00591    CubitBoolean is_pnt_on_ln( double pnt[3], double ln_orig[3], double ln_vec[3] );
00592
00596    CubitBoolean is_pnt_on_ln_seg( double pnt1[3], double end1[3], double end2[3] );
00597
00600    CubitBoolean is_pnt_on_pln( double pnt[3], double pln_orig[3], double pln_norm[3] );
00601
00605    CubitBoolean is_ln_on_pln( double ln_orig[3], double ln_vec[3],
00606                               double pln_orig[3], double pln_norm[3] );
00607
00608
00610    CubitBoolean is_pln_on_pln( double pln_orig1[3], double pln_norm1[3],
00611                                double pln_orig2[3], double pln_norm2[3] );
00612
00613    //**************************************************************************
00614    // Arcs/Circles
00615    //**************************************************************************
00616
00618
00619
00620
00621
00622
00623    void setup_arc( double vec_1[3], double vec_2[3], double origin[3],
00624                   double start_angle, double end_angle, // Angles in radians
00625                   double radius, AGT_3D_Arc &arc );
00626    void setup_arc( CubitVector& vec_1, CubitVector& vec_2,
00627                    CubitVector origin, double start_angle, // Angles in radians
00629                    AGT_3D_Arc &arc );
00631
00633
00634
00635    void get_arc_xyz( AGT_3D_Arc &arc, double param, double pnt[3] );
00636    void get_arc_xyz( AGT_3D_Arc &arc, double param, CubitVector& pnt );
00638
00645    int get_num_circle_tess_pnts( double radius, double len_tolerance = 1e-4 );
00646
00647    //**************************************************************************
00648    // Miscellaneous
00649    //**************************************************************************
00650
00653    void get_pln_orig_norm( double A, double B, double C, double D,
00654                            double pln_orig[3], double pln_norm[3] = NULL );
00658    void get_box_corners( double box_min[3],double box_max[3],double c[8][3] );
00659
00661
00662
00663
00664
00665
00666
00667
00668
00669
00670
00671    CubitStatus min_pln_box_int_corners( const CubitPlane& plane,
00672                                         const CubitBox& box,
00673                                         int extension_type, double extension,
00674                                         CubitVector& p1, CubitVector& p2,
00675                                         CubitVector& p3, CubitVector& p4,
00676                                         CubitBoolean silent = CUBIT_FALSE );
00677    CubitStatus min_pln_box_int_corners( CubitVector& vec1,
00678                                         CubitVector& vec2,
00679                                         CubitVector& vec3,
00680                                         CubitVector& box_min, CubitVector& box_max,
00681                                         int extension_type, double extension,
00682                                         CubitVector& p1, CubitVector& p2,
00683                                         CubitVector& p3, CubitVector& p4,
00684                                         CubitBoolean silent = CUBIT_FALSE );
00685    CubitStatus min_pln_box_int_corners( double plane_norm[3], double plane_coeff,
00686                                        double box_min[3], double box_max[3],
00687                                        int extension_type, double extension,
00688                                        double p1[3], double p2[3],   // OUT
00689                                        double p3[3], double p4[3],   // OUT
00690                                        CubitBoolean silent = CUBIT_FALSE);
00692
00694
00695
00696    CubitStatus get_tight_bounding_box( DLIList<CubitVector*> &point_list,
00697                                        CubitVector& p1, CubitVector& p2,
00698                                        CubitVector& p3, CubitVector& p4,
00699                                        CubitVector& p5, CubitVector& p6,
00700                                        CubitVector& p7, CubitVector& p8);
00701    CubitStatus get_tight_bounding_box( DLIList<CubitVector*> &point_list,
00702                                        CubitVector &center,
00703                                        CubitVector axes[3],
00704                                        CubitVector &extension );
00706
00707
00709
00710
00711
00712
00713
00714
00715
00716
00717
00718
00719    CubitStatus min_cyl_box_int( double radius, CubitVector& axis, CubitVector& center,
00720                                 CubitBox& box,
00721                                 int extension_type, double extension,
00722                                 CubitVector &start, CubitVector &end,
00723                                 int num_tess_pnts = 2048 );
00724    CubitStatus min_cyl_box_int( double radius, double axis[3], double center[3],
00725                                 double box_min[3], double box_max[3],
00726                                 int extension_type, double extension,
00727                                 double start[3], double end[3],
00728                                 int num_tess_pnts = 2048 );
00730
00732    double MinTriangleTriangle (Triangle3& tri0,
00733                                Triangle3& tri1,
00734                               double& s, double& t, double& u, double& v);
00735
00737    double AreaIntersection (Triangle &tri1, Triangle &tri2 );
00738
00740    double Angle( Triangle3& tri1, Triangle3& tri2 );
00741
00743    double MinPointTriangle (const Point3& p, const Triangle3& tri, double& s, double& t);
00744
00746    void Normal( Triangle3& tri, double normal[3] );
00747
00749    double ProjectedOverlap( Triangle3& tri1, Triangle3& tri2, bool draw_overlap = false );
00750
00751 protected:
00752    AnalyticGeometryTool();
00753    //- Class Constructor. (Not callable by user code. Class is constructed
00754    //- by the {instance()} member function.
00755
00756 private:
00757
00758    static AnalyticGeometryTool* instance_;
00759
00760    double agtEpsilon;  // default = 1e-6
00761
00762
00764 //                               MAGIC SOFTWARE
00765 //This code was obtained from Dave Eberly, at www.magic-software.com.  It
00766 //has been modified only to be placed into AnalyticGeometryTool.  This was
00767 //done for convenience only.  An alternate solution would be to create a
00768 //separate library for these functions.
00769 //
00770 //Steve Storm, storm@CAT.com, 3-27-99
00772 //MAGIC is an acronym for My Alternate Graphics and Image Code. The
00773 //initial code base originated in 1991 as an attempt to answer questions that
00774 //arise in my favorite news group, comp.graphics.algorithms, and has been
00775 //growing ever since. Magic is intended to provide free source code for solving
00776 //problems that commonly arise in computer graphics and image analysis. While
00777 //the code at this web site is free, additional conditions are:
00778 //
00779 // *  You may distribute the original source code to others at no charge. You
00780 //    got it for free, so do not charge others for it.
00781 // *  You may modify the original source code and distribute it to others at no
00782 //    charge. The modified code must be documented to indicate that it is not
00783 //    part of the original package. I do not want folks to blame me for bugs
00784 //    introduced by modifications. I do accept blame for bugs that are my
00785 //    doing and will do my best to fix them.
00786 // *  You may use this code for non-commercial purposes. You may also
00787 //    incorporate this code into commercial packages. However, you may not
00788 //    sell any of your source code which contains my original and/or modified
00789 //    source code (see items 1 and 2 in this list of conditions). In such a case,
00790 //    you would need to factor out my code and freely distribute it. Send me
00791 //    email if you use it. I am always interested in knowing where and how
00792 //    my code is being used.
00793 // *  The original code comes with absolutely no warranty and no guarantee is
00794 //    made that the code is bug-free. Caveat emptor.
00795 //
00796 // Dave Eberly, eberly@magic-software.com, www.magic-software.com
00798 //FILE: minbox2.h
00799    double Area (int N, Point2* pt, double angle);
00800    void MinimalBoxForAngle (int N, Point2* pt, double angle, OBBox2& box);
00801    OBBox2 MinimalBox2 (int N, Point2* pt);
00802    // Functions to find minimal area rectangle that surrounds a set of points.
00803
00804 #if 1
00805 //FILE: minbox3.h
00806    void MatrixToAngleAxis( double** R, double& angle, double axis[3] );
00807    void AngleAxisToMatrix( double angle, double axis[3], double R[3][3] );
00808    double Volume( int N, Point3* pt, double angle[3] );
00809    void MinimalBoxForAngles( int N, Point3* pt, double angle[3], OBBox3& box );
00810    void GetInterval( double A[3], double D[3], double& tmin, double& tmax );
00811    void Combine( double result[3], double A[3], double t, double D[3] );
00812    double MinimizeOnInterval( int N, Point3* pt, double A[3], double D[3] );
00813    double MinimizeOnLattice( int N, Point3* pt, double A[3], int layers,
00814                              double thickness );
00815    void InitialGuess( int N, Point3* pt, double angle[3] );
00816    OBBox3 MinimalBox3( int N, Point3* pt );
00817    // Functions to find minimal volume box that surrounds a set of points.
00818 #endif
00819
00820    double Abs (double x);
00821    double ACos (double x);
00822    double ATan2 (double y, double x);
00823    double Cos (double x);
00824    double Sign (double x);
00825    double Sin (double x);
00826    double Sqrt (double x);
00827    double UnitRandom ();
00828    double Dot (const Point2& p, const Point2& q);
00829    double Length (const Point2& p);
00830    CubitBoolean Unitize (Point2& p );
00831    double Dot (const Point3& p, const Point3& q);
00832    Point3 Cross (const Point3& p, const Point3& q);
00833    double Length (const Point3& p);
00834    CubitBoolean Unitize (Point3& p );
00835    // Supporting code for triangle calculations
00836
00837    double MinLineSegmentLineSegment (const Line3& seg0, const Line3& seg1, double& s, double& t);
00838
00839    //get intersection between bounding box and plane
00840    int get_plane_bbox_intersections( double box_min[3],
00841                                      double box_max[3],
00842                                      double pln_orig[3],
00843                                      double pln_norm[3],
00844                                      double int_array[6][3] );
00845 // FILE: pt3tri3.cpp
00846
00847
00848 // FILE: lin3tri3.cpp
00849    double MinLineSegmentTriangle (const Line3& seg, const Triangle3& tri, double& r, double& s, double& t);
00850    // Code to support distance between triangles
00851
00852 //FILE: triasect.cpp
00853    AgtLine EdgeToLine (Point2* v0, Point2* v1);
00854    void TriangleLines (Triangle* T, AgtLine line[3]);
00855    AgtTriList* SplitAndDecompose (Triangle* T, AgtLine* line);
00856    AgtTriList* Intersection (Triangle* T0, Triangle* T1);
00857    double AreaTriangle (Triangle* T);
00858    double Area (AgtTriList* list);
00859    double Area (Triangle &tri1, Triangle &tri2 );
00860    // Code to support area of overlap of two triangles
00861 };
00862
00863 #if 1
00864 // FILE: eigen.h
00865 class mgcEigen
00866 {
00867 public:
00868         mgcEigen (int _size);
00869         ~mgcEigen ();
00870
00871         mgcEigen& Matrix (float** inmat);
00872
00873 private:
00874         int size;
00875         float** mat;
00876         float* diag;
00877         float* subd;
00878
00879         // Householder reduction to tridiagonal form
00880         void Tridiagonal2 (float** pmat, float* pdiag, float* psubd);
00881         void Tridiagonal3 (float** pmat, float* pdiag, float* psubd);
00882         void Tridiagonal4 (float** pmat, float* pdiag, float* psubd);
00883         void TridiagonalN (int n, float** mat, float* diag, float* subd);
00884
00885         // QL algorithm with implicit shifting, applies to tridiagonal matrices
00886         void QLAlgorithm (int n, float* pdiag, float* psubd, float** pmat);
00887
00888         // sort eigenvalues from largest to smallest
00889         void DecreasingSort (int n, float* eigval, float** eigvec);
00890
00891         // sort eigenvalues from smallest to largest
00892         void IncreasingSort (int n, float* eigval, float** eigvec);
00893
00894 // error handling
00895 public:
00896         static int verbose1;
00897         static unsigned error;
00898         static void Report (std::ostream& ostr);
00899 private:
00900         static const unsigned invalid_size;
00901         static const unsigned allocation_failed;
00902         static const unsigned ql_exceeded;
00903         static const char* message[3];
00904         static int Number (unsigned single_error);
00905         static void Report (unsigned single_error);
00906 };
00907
00908 // FILE: eigen.h
00909 class mgcEigenD
00910 {
00911 public:
00912         mgcEigenD (int _size);
00913         ~mgcEigenD ();
00914
00915     // set the matrix for eigensolving
00916         double& Matrix (int row, int col) { return mat[row][col]; }
00917         mgcEigenD& Matrix (double** inmat);
00918
00919         // get the results of eigensolving
00920         double Eigenvector (int row, int col) { return mat[row][col]; }
00921         const double** Eigenvector () { return (const double**) mat; }
00922
00923         void EigenStuff3 ();  // uses TriDiagonal3
00924
00925 private:
00926         int size;
00927         double** mat;
00928         double* diag;
00929         double* subd;
00930
00931         // Householder reduction to tridiagonal form
00932         void Tridiagonal2 (double** mat, double* diag, double* subd);
00933         void Tridiagonal3 (double** mat, double* diag, double* subd);
00934         void Tridiagonal4 (double** mat, double* diag, double* subd);
00935         void TridiagonalN (int n, double** mat, double* diag, double* subd);
00936
00937         // QL algorithm with implicit shifting, applies to tridiagonal matrices
00938         void QLAlgorithm (int n, double* diag, double* subd, double** mat);
00939
00940         // sort eigenvalues from largest to smallest
00941         void DecreasingSort (int n, double* eigval, double** eigvec);
00942
00943         // sort eigenvalues from smallest to largest
00944         void IncreasingSort (int n, double* eigval, double** eigvec);
00945
00946 // error handling
00947 public:
00948         static int verbose1;
00949         static unsigned error;
00950         static void Report (std::ostream& ostr);
00951 private:
00952         static const unsigned invalid_size;
00953         static const unsigned allocation_failed;
00954         static const unsigned ql_exceeded;
00955         static const char* message[3];
00956         static int Number (unsigned single_error);
00957         static void Report (unsigned single_error);
00958 };
00959 #endif
00960
00962 //---------------------------------------------------------------------------
00963 // Wrapped math functions
00964 //---------------------------------------------------------------------------
00965 inline double AnalyticGeometryTool::Abs (double x)
00966 {
00967     return double(fabs(x));
00968 }
00969
00970 inline double AnalyticGeometryTool::ATan2 (double y, double x)
00971 {
00972     return double(atan2(y,x));
00973 }
00974
00975 inline double AnalyticGeometryTool::Sqrt (double x)
00976 {
00977     return double(sqrt(x));
00978 }
00979
00980 inline double AnalyticGeometryTool::UnitRandom ()
00981 {
00982     return double(rand())/double(RAND_MAX);
00983 }
00984
00985 //---------------------------------------------------------------------------
00986 // Points in 2D
00987 //---------------------------------------------------------------------------
00988
00989 inline double AnalyticGeometryTool::Dot (const Point3& p, const Point3& q)
00990 {
00991     return double(p.x*q.x + p.y*q.y + p.z*q.z);
00992 }
00993
00994
00995 inline double AnalyticGeometryTool::set_epsilon( double new_epsilon )
00996 {
00997    double old_epsilon = agtEpsilon;
00998    agtEpsilon = new_epsilon;
00999    return old_epsilon;
01000 }
01001
01002 inline CubitBoolean AnalyticGeometryTool::dbl_eq(double val_1,double val_2)
01003 {
01004    CubitBoolean result;
01005    if (fabs(val_1 - val_2) < agtEpsilon)
01006       result = CUBIT_TRUE;
01007    else
01008       result = CUBIT_FALSE;
01009    return result;
01010 }
01011
01012 #endif
01013
```