Main Page | Namespace List | Class Hierarchy | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

Geometry.h

Go to the documentation of this file.
00001 /*
00002     liban - base function to do usefull things related to simulation and scientific work
00003     Copyright (C) 1999-2005 Stephane Magnenat <stephane at magnenat dot net>
00004     Copyright (c) 2004-2005 Antoine Beyeler <antoine dot beyeler at epfl dot ch>
00005     Copyright (C) 2005 Laboratory of Intelligent Systems, EPFL, Lausanne
00006     See AUTHORS for details
00007 
00008     This program is free software. You can redistribute it and/or modify
00009     it under the terms of the GNU General Public License as published by
00010     the Free Software Foundation; either version 2 of the License, or
00011     (at your option) any later version.
00012 
00013     This program is distributed in the hope that it will be useful,
00014     but WITHOUT ANY WARRANTY; without even the implied warranty of
00015     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00016     GNU General Public License for more details.
00017 
00018     You should have received a copy of the GNU General Public License
00019     along with this program; if not, write to the Free Software
00020     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00021 */
00022 
00023 #ifndef __AN_GEOMETRY_H
00024 #define __AN_GEOMETRY_H
00025 
00026 #ifdef WIN32
00027 #define _USE_MATH_DEFINES
00028 #endif
00029 #include <cmath>
00030 #include <vector>
00031 
00036 
00037 namespace An
00038 {
00040 
00047     struct Vector
00048     {
00050         double x, y; 
00051     
00053         Vector() { x = y = 0; }
00055         Vector(double v) { this->x = v; this->y = v; }
00057         Vector(double x, double y) { this->x = x; this->y = y; }
00059         Vector(double array[2]) { x = array[0]; y = array[1]; }
00060     
00062         void operator +=(const Vector &v) { x += v.x; y += v.y; }
00064         void operator -=(const Vector &v) { x -= v.x; y -= v.y; }
00066         void operator *=(double f) { x *= f; y *= f; }
00068         void operator /=(double f) { x /= f; y /= f; }
00070         Vector operator +(const Vector &v) const { Vector n; n.x = x + v.x; n.y = y + v.y; return n; }
00072         Vector operator -(const Vector &v) const { Vector n; n.x = x - v.x; n.y = y - v.y; return n; }
00074         Vector operator /(double f) const { Vector n; n.x = x/f; n.y = y/f; return n; }
00076         Vector operator *(double f) const { Vector n; n.x = x*f; n.y = y*f; return n; }
00077     
00079         double operator *(const Vector &v) const { return x*v.x + y*v.y; }
00081         double norm(void) const { return sqrt(x*x + y*y); }
00083         double norm2(void) const { return x*x+y*y; }
00085         double cross(const Vector &v) const { return x * v.y - y * v.x; }
00087         Vector unitary(void) const { Vector n; if (norm()) { n = *this; n /= norm(); } return n; }
00089         double angle(void) const { return atan2(y, x); }
00090     };
00091     
00093 
00094     typedef Vector Point;
00095     
00097 
00104     struct Matrix22
00105     {
00107         double a, b, c, d;
00108     
00110         Matrix22() { a = b = c = d = 0; }
00112         Matrix22(double a, double b, double c, double d) { this->a = a; this->b = b; this->c = c; this->d = d; }
00114         Matrix22(double alpha) { a = cos(alpha); b = sin(alpha); c = -b; d = a; }
00116         Matrix22(double array[4]) { a=array[0]; b=array[1]; c=array[2]; d=array[3]; }
00117     
00119         void operator +=(const Matrix22 &v) { a += v.a; b += v.b; c += v.c; d += v.d; }
00121         void operator -=(const Matrix22 &v) { a -= v.a; b -= v.b; c -= v.c; d -= v.d; }
00123         void operator *=(double f) { a *= f; b *= f; c *= f; d *= f; }
00125         void operator /=(double f) { a /= f; b /= f; c /= f; d /= f; }
00127         Matrix22 operator +(const Matrix22 &v) const { Matrix22 n; n.a = a + v.a; n.b = b + v.b; n.c = c + v.c; n.d = d + v.d; return n; }
00129         Matrix22 operator -(const Matrix22 &v) const { Matrix22 n; n.a = a - v.a; n.b = b - v.b; n.c = c - v.c; n.d = d - v.d; return n; }
00131         Matrix22 operator *(double f) const { Matrix22 n; n.a = a * f; n.b = b * f; n.c = c * f; n.d = d * f; return n; }
00133         Matrix22 operator /(double f) const { Matrix22 n; n.a = a / f; n.b = b / f; n.c = c / f; n.d = d / f; return n; }
00134         
00136         Point operator*(const Point &v) { Point n; n.x = v.x*a + v.y*c; n.y = v.x*b + v.y*d; return n; }
00137     };
00138     
00140 
00141     struct Segment
00142     {
00144         Segment(double ax, double ay, double bx, double by) { this->a.x = ax; this->a.y = ay; this->b.x = bx; this->b.y = by; }
00146         Segment(double array[4]) { a.x = array[0]; a.y = array[1]; b.x = array[2]; b.y = array[3]; }
00148         Segment(const Point &p1, const Point &p2) { a = p1; b = p2; }
00149         
00151         Point a, b;
00152     
00154         double dist(const Point &p) const
00155         {
00156             Vector n(a.y-b.y, b.x-a.x);
00157             Vector u = n.unitary();
00158             Vector ap = p-a;
00159             return ap * u;
00160         }
00161     
00163         bool doesIntersect(const Segment &o) const
00164         {
00165             double s2da = dist (o.a);
00166             double s2db = dist (o.b);
00167             double s1da = o.dist (a);
00168             double s1db = o.dist (b);
00169             return (s2da*s2db<0) && (s1da*s1db<0);
00170         }
00171     };
00172     
00174 
00175     typedef std::vector<Point> Polygone;
00176     
00178 
00179     inline double normalizeAngle(double angle)
00180     {
00181         while (angle > M_PI)
00182             angle -= 2*M_PI;
00183         while (angle < -M_PI)
00184             angle += 2*M_PI;
00185         return angle;
00186     }
00187   
00191 
00192     inline Point getIntersection(const Segment &s1, const Segment &s2)
00193     {
00194         // compute first segment's equation
00195         double c1 = s1.a.y + (-s1.a.x / (s1.b.x - s1.a.x)) * (s1.b.y - s1.a.y);
00196         double m1 = (s1.b.y - s1.a.y) / (s1.b.x - s1.a.x);
00197  
00198         // compute second segment's equation
00199         double c2 = s2.a.y + (-s2.a.x / (s2.b.x - s2.a.x)) * (s2.b.y - s2.a.y);
00200         double m2 = (s2.b.y - s2.a.y) / (s2.b.x - s2.a.x);
00201 
00202         // are the lines parallel ?
00203         if (m1 == m2)
00204             return Point(HUGE_VAL, HUGE_VAL);
00205 
00206         double x1 = s1.a.x;
00207         double x2 = s1.b.x;
00208         double x3 = s2.a.x;
00209         double x4 = s2.b.x;
00210         double y1 = s1.a.y;
00211         double y2 = s1.b.y;
00212         double y3 = s2.a.y;
00213         double y4 = s2.b.y;
00214 
00215         // make sure x1 < x2
00216         if (x1 > x2)
00217         {
00218             double temp = x1;
00219             x1 = x2;
00220             x2 = temp;
00221         }
00222 
00223         // make sure x3 < x4
00224         if (x3 > x4)
00225         {
00226             double temp = x3;
00227             x3 = x4;
00228             x4 = temp;
00229         }
00230 
00231         // make sure y1 < y2
00232         if (y1 > y2)
00233         {
00234             double temp = y1;
00235             y1 = y2;
00236             y2 = temp;
00237         }
00238 
00239         // make sure y3 < y4
00240         if (y3 > y4)
00241         {
00242             double temp = y3;
00243             y3 = y4;
00244             y4 = temp;
00245         }
00246 
00247         // intersection point in case of infinite slopes
00248         double x;
00249         double y;
00250 
00251         // infinite slope m1
00252         if (x1 == x2)
00253         {
00254             x = x1;
00255             y = m2 * x1 + c2;
00256             if (x > x3 && x < x4 && y > y1 && y <y2)
00257                 return Point(x, y);
00258             else
00259                 return Point(HUGE_VAL, HUGE_VAL);
00260         }
00261 
00262         // infinite slope m2
00263         if (x3 == x4)
00264         {
00265             x = x3;
00266             y = m1 * x3 + c1;
00267             if (x > x1 && x < x2 && y > y3 && y < y4)
00268                 return Point(x, y);
00269             else
00270                 return Point(HUGE_VAL, HUGE_VAL);
00271         }
00272         
00273         // compute lines intersection point
00274         x = (c2 - c1) / (m1 - m2);
00275 
00276         // see whether x in in both ranges [x1, x2] and [x3, x4]
00277         if (x > x1 && x < x2 && x > x3 && x < x4)
00278             return Point(x, m1 * x + c1);
00279   
00280         return Point(HUGE_VAL, HUGE_VAL);
00281     }
00282 }
00283 
00284 #endif

Generated on Mon Oct 24 17:30:33 2005 for liban by  doxygen 1.4.2