Back to index

texmacs  1.0.7.15
point.cpp
Go to the documentation of this file.
00001 
00002 /******************************************************************************
00003 * MODULE     : point.cpp
00004 * DESCRIPTION: points
00005 * COPYRIGHT  : (C) 2003  Joris van der Hoeven
00006 *******************************************************************************
00007 * This software falls under the GNU general public license version 3 or later.
00008 * It comes WITHOUT ANY WARRANTY WHATSOEVER. For details, see the file LICENSE
00009 * in the root directory or <http://www.gnu.org/licenses/gpl-3.0.html>.
00010 ******************************************************************************/
00011 
00012 #include "point.hpp"
00013 #include "math_util.hpp"
00014 
00015 point
00016 operator - (point p) {
00017   int i, n= N(p);
00018   point r (n);
00019   for (i=0; i<n; i++)
00020     r[i]= - p[i];
00021   return r;
00022 }
00023 
00024 point
00025 operator + (point p1, point p2) {
00026   int i, n= min (N(p1), N(p2));
00027   point r (n);
00028   for (i=0; i<n; i++)
00029     r[i]= p1[i] + p2[i];
00030   return r;
00031 }
00032 
00033 point
00034 operator - (point p1, point p2) {
00035   int i, n= min (N(p1), N(p2));
00036   point r (n);
00037   for (i=0; i<n; i++)
00038     r[i]= p1[i] - p2[i];
00039   return r;
00040 }
00041 
00042 point
00043 operator * (double x, point p) {
00044   int i, n= N(p);
00045   point r (n);
00046   for (i=0; i<n; i++)
00047     r[i]= x * p[i];
00048   return r;
00049 }
00050 
00051 point
00052 operator / (point p, double x) {
00053   int i, n= N(p);
00054   point r (n);
00055   for (i=0; i<n; i++)
00056     r[i]= p[i] / x;
00057   return r;
00058 }
00059 
00060 bool
00061 operator == (point p1, point p2) {
00062   if (N(p1) != N(p2)) return false;
00063   int i, n= N(p1);
00064   for (i=0; i<n; i++)
00065     if (!fnull (p1[i]-p2[i], 1e-6)) return false;
00066   return true;
00067 }
00068 
00069 bool
00070 is_point (tree t) {
00071   return L(t)==_POINT;
00072 }
00073 
00074 point
00075 as_point (tree t) {
00076   if (!is_tuple (t) && !is_point (t)) return point ();
00077   else {
00078     int i, n= N(t);
00079     point p(n);
00080     for (i=0; i<n; i++)
00081       p[i]= as_double (t[i]);
00082     return p;
00083   }
00084 }
00085 
00086 tree
00087 as_tree (point p) {
00088   int i, n= N(p);
00089   tree t (_POINT, n);
00090   for (i=0; i<n; i++)
00091     t[i]= as_string (p[i]);
00092   return t;
00093 }
00094 
00095 double
00096 operator * (point p1, point p2) {
00097   int i, n= min (N(p1), N(p2));
00098   double r= 0;
00099   for (i=0; i<n; i++)
00100     r+= p1[i] * p2[i];
00101   return r;
00102 }
00103 
00104 static point
00105 mult (double re, double im, point p) {
00106   if (N(p)==0) p= point (0.0, 0.0);
00107   if (N(p)==1) p= point (p[0], 0.0);
00108   return point (re * p[0] - im * p[1],
00109               re * p[1] + im * p[0]);
00110 }
00111 
00112 point
00113 rotate_2D (point p, point o, double angle) {
00114   return mult (cos (angle), sin (angle), p - o) + o;
00115 }
00116 
00117 double
00118 norm (point p) {
00119   return sqrt (p*p);
00120 }
00121 
00122 double
00123 arg (point p) {
00124   double n= norm(p);
00125   p=p/n;
00126   if (p[1]<0) return 2*tm_PI-acos(p[0]);
00127   else return acos(p[0]);
00128 }
00129 
00130 point
00131 proj (axis ax, point p) {
00132   int i, n= min (N(ax.p0), N(ax.p1));
00133   point a (n), b (n);
00134   for (i=0; i<n ; i++) {
00135     a[i]= ax.p1[i] - ax.p0[i];
00136     b[i]= ax.p0[i];
00137   }
00138   if (norm (a) < 1.0e-6)
00139     return ax.p0;
00140   else
00141     return b + ((a*p - a*b) / (a*a)) * a;
00142 }
00143 
00144 double
00145 dist (axis ax, point p) {
00146   return norm (p - proj (ax, p));
00147 }
00148 
00149 double
00150 seg_dist (axis ax, point p) {
00151   point ab= ax.p1 - ax.p0;
00152   point ba= ax.p0 - ax.p1;
00153   point ap= p - ax.p0;
00154   point bp= p - ax.p1;
00155   if (ab * ap > 0 && ba * bp > 0)
00156     return dist (ax, p);
00157   else
00158     return min (norm (ap), norm (bp));
00159 }
00160 
00161 bool
00162 collinear (point p1, point p2) {
00163   return fnull (fabs (p1*p2) - norm(p1)*norm(p2), 1.0e-6);
00164 }
00165 
00166 bool
00167 linearly_dependent (point p1, point p2, point p3) {
00168   return fnull (norm (p1-p2), 1e-6) ||
00169         fnull (norm (p2-p3), 1e-6) ||
00170         fnull (norm (p3-p1), 1e-6) ||
00171         collinear (p2-p1, p3-p1);
00172 }
00173 
00174 bool orthogonalize (point &i, point &j, point p1, point p2, point p3) {
00175   if (linearly_dependent (p1, p2, p3)) return false;
00176   i= (p2-p1) / norm (p2-p1);
00177   j= (p3-p1) - ((p3-p1) * i) * i;
00178   j= j / norm (j);
00179   return true;
00180 }
00181 
00182 axis
00183 midperp (point p1, point p2, point p3) {
00184   axis a;
00185   if (linearly_dependent (p1, p2, p3))
00186     a.p0= a.p1= point (0);
00187   else {
00188     point i, j;
00189     orthogonalize (i, j, p1, p2, p3);
00190     a.p0= (p1+p2) / 2;
00191     a.p1= a.p0 + j;
00192   }
00193   return a;
00194 }
00195 
00196 point
00197 intersection (axis A, axis B) {
00198   point i, j;
00199   if (!orthogonalize (i, j, A.p0, A.p1, B.p0)) {
00200     if (orthogonalize (i, j, A.p0, A.p1, B.p1))
00201       return B.p0;
00202     else
00203       return point (0);
00204   }
00205   point a(2), b(2), u(2), v(2), p(2);
00206   a[0]= a[1]= 0;
00207   u[0]= (A.p1 - A.p0) * i;
00208   u[1]= (A.p1 - A.p0) * j;
00209   b[0]= (B.p0 - A.p0) * i;
00210   b[1]= (B.p0 - A.p0) * j;
00211   v[0]= (B.p1 - B.p0) * i;
00212   v[1]= (B.p1 - B.p0) * j;
00213   if (fnull (norm (u), 1e-6) ||
00214       fnull (norm (v), 1e-6) ||
00215       collinear (u, v))
00216     return point (0);
00217   else {
00218     double t;
00219     t= (v[0] * (b[1]-a[1]) + v[1] * (a[0]-b[0]))
00220        /
00221        (v[0]*u[1] - v[1]*u[0]);
00222     return A.p0 + t * (u[0]*i + u[1]*j);
00223   }
00224 }
00225 
00226 bool
00227 inside_rectangle (point p, point p1, point p2) {
00228   return p[0]>=p1[0] && p[1]>=p1[1]
00229       && p[0]<=p2[0] && p[1]<=p2[1];
00230 }