src/geom.cpp

   1 /**
   2  *  \file src/geom.cpp
   3  *  \brief Various geometrical calculations.
   4  */
   5
   6 #ifdef HAVE_CONFIG_H
   7 # include <config.h>
   8 #endif
   9 #include "geom.h"
  10 #include <libnr/nr-point-fns.h>
  11
  12 /**
  13  * Finds the intersection of the two (infinite) lines
  14  * defined by the points p such that dot(n0, p) == d0 and dot(n1, p) == d1.
  15  *
  16  * If the two lines intersect, then \a result becomes their point of
  17  * intersection; otherwise, \a result remains unchanged.
  18  *
  19  * This function finds the intersection of the two lines (infinite)
  20  * defined by n0.X = d0 and x1.X = d1.  The algorithm is as follows:
  21  * To compute the intersection point use Cramer's rule:
  22  * (see http://en.wikipedia.org/wiki/Cramer%27s_rule)
  23  * \verbatim
  24  * convert lines to form
  25  * ax + by = c
  26  * dx + ey = f
  27  *
  28  * (
  29  *  e.g. a = (x2 - x1), b = (y2 - y1), c = (x2 - x1)*x1 + (y2 - y1)*y1
  30  * )
  31  *
  32  * In our case we use:
  33  *   a = n0.x     d = n1.x
  34  *   b = n0.y     e = n1.y
  35  *   c = d0        f = d1
  36  *
  37  * so:
  38  *
  39  * adx + bdy = cd
  40  * adx + aey = af
  41  *
  42  * bdy - aey = cd - af
  43  * (bd - ae)y = cd - af
  44  *
  45  * y = (cd - af)/(bd - ae)
  46  *
  47  * repeat for x and you get:
  48  *
  49  * x = (fb - ce)/(bd - ae)                \endverbatim
  50  *
  51  * If the denominator (bd-ae) is 0 then the lines are parallel, if the
  52  * numerators are then 0 then the lines coincide.
  53  *
  54  * \todo Why not use existing but outcommented code below
  55  * (HAVE_NEW_INTERSECTOR_CODE)?
  56  */
  57 IntersectorKind intersector_line_intersection(NR::Point const &n0, double const d0,
  58                                               NR::Point const &n1, double const d1,
  59                                               NR::Point &result)
  60 {
  61     double denominator = dot(NR::rot90(n0), n1);
  62     double X = n1[NR::Y] * d0 -
  63                n0[NR::Y] * d1;
  64     /* X = (-d1, d0) dot (n0[Y], n1[Y]) */
  65
  66     if (denominator == 0) {
  67         if ( X == 0 ) {
  68             return COINCIDENT;
  69         } else {
  70             return PARALLEL;
  71         }
  72     }
  73
  74     double Y = n0[NR::X] * d1 -
  75                n1[NR::X] * d0;
  76
  77     result = NR::Point(X, Y) / denominator;
  78
  79     return INTERSECTS;
  80 }
  81
  82
  83
  84
  85 /*
  86  New code which we are not yet using
  87 */
  88 #ifdef HAVE_NEW_INTERSECTOR_CODE
  89
  90
  91 /* ccw exists as a building block */
  92 static int
  93 sp_intersector_ccw(const NR::Point p0, const NR::Point p1, const NR::Point p2)
  94 /* Determine which way a set of three points winds. */
  95 {
  96         NR::Point d1 = p1 - p0;
  97         NR::Point d2 = p2 - p0;
  98 /* compare slopes but avoid division operation */
  99         double c = dot(NR::rot90(d1), d2);
 100         if(c > 0)
 101                 return +1; // ccw - do these match def'n in header?
 102         if(c < 0)
 103                 return -1; // cw
 104
 105         /* Colinear [or NaN].  Decide the order. */
 106         if ( ( d1[0] * d2[0] < 0 )  ||
 107              ( d1[1] * d2[1] < 0 ) ) {
 108                 return -1; // p2  <  p0 < p1
 109         } else if ( dot(d1,d1) < dot(d2,d2) ) {
 110                 return +1; // p0 <= p1  <  p2
 111         } else {
 112                 return 0; // p0 <= p2 <= p1
 113         }
 114 }
 115
 116 /** Determine whether two line segments intersect.  This doesn't find
 117     the point of intersection, use the line_intersect function above,
 118     or the segment_intersection interface below.
 119
 120     \pre neither segment is zero-length; i.e. p00 != p01 and p10 != p11.
 121  */
 122 static bool
 123 sp_intersector_segment_intersectp(NR::Point const &p00, NR::Point const &p01,
 124                                   NR::Point const &p10, NR::Point const &p11)
 125 {
 126         g_return_val_if_fail(p00 != p01, false);
 127         g_return_val_if_fail(p10 != p11, false);
 128
 129         /* true iff (    (the p1 segment straddles the p0 infinite line)
 130          *           and (the p0 segment straddles the p1 infinite line) ). */
 131         return ((sp_intersector_ccw(p00,p01, p10)
 132                  *sp_intersector_ccw(p00, p01, p11)) <=0 )
 133                 &&
 134                 ((sp_intersector_ccw(p10,p11, p00)
 135                   *sp_intersector_ccw(p10, p11, p01)) <=0 );
 136 }
 137
 138
 139 /** Determine whether \& where two line segments intersect.
 140
 141     If the two segments don't intersect, then \a result remains unchanged.
 142
 143     \pre neither segment is zero-length; i.e. p00 != p01 and p10 != p11.
 144 **/
 145 static sp_intersector_kind
 146 sp_intersector_segment_intersect(NR::Point const &p00, NR::Point const &p01,
 147                                  NR::Point const &p10, NR::Point const &p11,
 148                                  NR::Point &result)
 149 {
 150         if(sp_intersector_segment_intersectp(p00, p01, p10, p11)) {
 151                 NR::Point n0 = (p00 - p01).ccw();
 152                 double d0 = dot(n0,p00);
 153
 154                 NR::Point n1 = (p10 - p11).ccw();
 155                 double d1 = dot(n1,p10);
 156                 return sp_intersector_line_intersection(n0, d0, n1, d1, result);
 157         } else {
 158                 return no_intersection;
 159         }
 160 }
 161
 162 #endif /* end yet-unused HAVE_NEW_INTERSECTOR_CODE code */
 163
 164 /*
 165   Local Variables:
 166   mode:c++
 167   c-file-style:"stroustrup"
 168   c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
 169   indent-tabs-mode:nil
 170   fill-column:99
 171   End:
 172 */
 173 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=99 :