home *** CD-ROM | disk | FTP | other *** search
/ vis-ftp.cs.umass.edu / vis-ftp.cs.umass.edu.tar / vis-ftp.cs.umass.edu / pub / Software / ASCENDER / ascendMar8.tar / UMass / BuildingFinder / Staging / intersect.c < prev    next >
C/C++ Source or Header  |  1995-05-09  |  5KB  |  277 lines
  1. /*
  2.  *
  3.  * IntersectLines
  4.  *
  5.  * Compute the Intersection of two lines using the Parametric 
  6.  * representations.
  7.  *
  8.  *
  9.  *    _                       _    _   _        _       _
  10.  *   | (x2 - x1)     (x3 - x4) |  |  t  |      | x3 - x1 | 
  11.  *   | (y2 - y1)     (y3 - y4) |  |  s  |  =   | y3 - y1 | 
  12.  *    -                       -    -   -        -       -
  13.  *
  14.  */
  15. #include "../polygons.h"
  16.  
  17. #define T        2
  18. #define L        3
  19. #define N        -1
  20. #define LOC_PERCENTAGE    0.8
  21. #define PIBY2  1.570796
  22. #define BIG_THETA 0.2
  23. #define MIN_ANGLE_THRESH 1.2
  24.  
  25.  
  26. #define SMALL(x,y)          ((x < y) ? x : y)
  27. #define nearEnd(x,val)          ((fabs(x) >= val) || (fabs(1.0 - x) >= val))
  28. #define onLine(t,val)               ((t >= (1.0 - val)) && (t <= val))
  29.  
  30.  
  31.  
  32.  
  33. int conflictIntersection(int x1, int y1, int x2, int y2)
  34. /*
  35.  * Given endpoints of two lines, is there an intersection that occurs 
  36.  * between the line and any other lines in the line list.
  37.  *
  38.  *
  39.  */
  40. {
  41.  
  42.     double_2 *lines_form_corner();
  43.     line *lptr;
  44.     float val;
  45.     line lrec;
  46.     double parms[2];
  47.  
  48.     double_2 Ipoint;
  49.     double_2 *apoint;
  50.  
  51.     lrec.x1 = (double) x1; lrec.y1 = (double) y1;
  52.     lrec.x2 = (double) x2; lrec.y2 = (double) y2;
  53.  
  54.  
  55.  
  56.  
  57.     lptr = line_list;
  58.     while (lptr != NULL) {
  59.  
  60.        intersect_point(lrec,*lptr, &Ipoint);
  61.        computeParms(Ipoint,lrec,*lptr,parms);
  62.  
  63.        if (onLine(parms[0],LOC_PERCENTAGE) && 
  64.         onLine(parms[1],LOC_PERCENTAGE) && angleBetween(lrec,*lptr) >
  65.         MIN_ANGLE_THRESH ) {
  66.           return(TRUE);
  67.             
  68.        }
  69.  
  70.        lptr = lptr->next;
  71.     }
  72.  
  73.     return(FALSE);
  74. }
  75.  
  76.  
  77.  
  78. /***
  79. float IntersectLines(line l1, line l2, double Parm[2])
  80. {
  81.  
  82.  
  83.     float **A;
  84.     float *B;
  85.     int *index;
  86.     float d;
  87.     int j;
  88.     Point p;
  89.  
  90.     A = matrix(1,2,1,2);
  91.     B = vector(1,2);
  92.     index = ivector(1,2);
  93.  
  94.  
  95.     fprintf(stderr,"Intersect lines: (%f %f, %f %f) -- (%f %f, %f %f)\n",
  96.             l1.x1,l1.y1,l1.x2,l1.y2,l2.x1,l2.y1,l2.x2,l2.y2);
  97.  
  98.     A[1][1] = (float)(l1.x2 - l1.x1);
  99.     A[1][2] = (float)(l2.x1 - l2.x2);
  100.     A[2][1] = (float)(l1.y2 - l1.y1);
  101.     A[2][2] = (float)(l2.y1 - l2.y2);
  102.     
  103.     B[1] = (l2.x1 - l1.x1);
  104.     B[2] = (l2.y1 - l1.y1);
  105.  
  106.  
  107.     ludcmp(A, 2, index, &d);
  108.     for (j=1; j <= 2; j++) 
  109.         d *= A[j][j];
  110.     if (d == 0) {
  111.         Parm[0] = -1.0;
  112.         Parm[1] = -1.0;
  113.         return(-1.0);
  114.     }
  115.     lubksb(A, 2, index, B);
  116.  
  117.  
  118.         Parm[0] = B[1];
  119.     Parm[1] = B[2];
  120.  
  121.     free_matrix(A,1.0,2.0,1.0,2.0);
  122.     free_vector(B,1.0,2.0);
  123.  
  124.     return(Parm[0]);
  125. }
  126. ***/
  127.  /*
  128.  
  129. double
  130. angleBetween(line l1, line l2)
  131.  * 
  132.  * Compute the angle between two lines as:
  133.  *
  134. {
  135.     double theta1, theta2;
  136.  
  137.     theta1 = LineSegTheta(l1.x1,l1.y1,l1.x2,l1.y2);
  138.     theta2 = LineSegTheta(l2.x1,l2.y1,l2.x2,l2.y2);
  139.  
  140.  
  141.     return(  fabs(theta1 - theta2));
  142. }
  143. */
  144.  
  145.  
  146. double lineTheta(double StartCol, double StartRow, double EndCol, double EndRow)
  147. {
  148.  
  149.  double dx,dy;             /* Buffers */
  150.  
  151.  double Orientation;
  152.  
  153.   dx = EndCol - StartCol;
  154.   dy = EndRow - StartRow;
  155.  
  156.   if (dx == 0) return(PIBY2);     /* Line is vertical */
  157.   if (dy == 0) return(0);         /* Line is Horizontal */
  158.  
  159.   Orientation =  atan2( dy, dx );
  160.  
  161. /*
  162.  Make it so line orientation is always between 0 and pi.
  163.  
  164. */
  165.  
  166.   if (Orientation < 0 )
  167.       Orientation = Orientation + PI;
  168.  
  169.   if (Orientation > PI)
  170.       Orientation = Orientation - PI;
  171.  
  172.   return(Orientation);
  173.  
  174. }
  175.  
  176. int IntersectType(double t, double s)
  177. {
  178.  
  179.  
  180.  
  181.  
  182.         if (nearEnd(t,LOC_PERCENTAGE) && nearEnd(s,LOC_PERCENTAGE))
  183.                 return(L);
  184.         if (nearEnd(t,LOC_PERCENTAGE) && onLine(s,LOC_PERCENTAGE))
  185.                 return(T);
  186.         if (onLine(t,LOC_PERCENTAGE) && nearEnd(s,LOC_PERCENTAGE))
  187.                 return(T);
  188.         if (onLine(t,LOC_PERCENTAGE) && onLine(s,LOC_PERCENTAGE))
  189.                 return(X);
  190.  
  191.         return(N);
  192. }
  193.  
  194.  
  195. /* 
  196.  *
  197.  * Intersect
  198.  *
  199.  */
  200. void
  201. intersect_point(line l1, line l2, double_2 *result)
  202. {
  203.  
  204.     Point vec1, vec2;
  205.     double scale;
  206.     double a1, b1, c1, a2, b2, c2;
  207.     double div;
  208.     
  209.     a1 = l1.y1 - l1.y2;
  210.     b1 = l1.x2 - l1.x1;
  211.     c1 = (l1.x1 * l1.y2) - (l1.x2 * l1.y1);
  212.     scale = sqrt ( SQ(a1) + SQ(b1) );
  213.     if (scale == 0)  {
  214.         result->d1 = -2.0;
  215.         result->d2 = -2.0;
  216.         return;
  217.     }
  218.     
  219.     vec1.x = a1 / scale; vec1.y = b1 /scale; vec1.z = c1 /scale;
  220.  
  221.     a2 = l2.y1 - l2.y2;
  222.     b2 = l2.x2 - l2.x1;
  223.     c2 = (l2.x1 * l2.y2) - (l2.x2 * l2.y1);
  224.     scale = sqrt ( SQ(a2) + SQ(b2) );
  225.     if (scale == 0.0) {
  226.         result->d1 = -2.0;
  227.         result->d2 = -2.0;
  228.         return;
  229.     }
  230.  
  231.     vec1.x = a2 / scale; vec1.y = b2 /scale; vec1.z = c2 /scale;
  232.  
  233.  
  234.     div = (a1 * b2) - (a2 * b1);
  235.     if (div == 0.0) {
  236.         result->d1 = -2.0;
  237.         result->d2 = -2.0;
  238.         return;
  239.     }
  240.  
  241.     result->d1 = ((b1 * c2) - (b2 * c1)) / div;
  242.     result->d2 = ((c1 * a2) - (c2 * a1)) / div;
  243. }
  244.  
  245.  
  246. void
  247. computeParms(double_2 vert, line l1, line l2, double parms[2])
  248. {
  249.  
  250.         Point vec1, vec2;
  251.  
  252.         vec1.x  = l1.x2 - l1.x1;
  253.         vec1.y = l1.y2 - l1.y1;
  254.  
  255.         vec2.x = vert.d1 - l1.x1;
  256.         vec2.y = vert.d2 - l1.y1;
  257.  
  258.         if (vec1.x == 0.0 && vec1.y == 0.0) /* Zero length !! */
  259.                 parms[0] = -1.0;
  260.         else
  261.                 parms[0] = (vec1.x*vec2.x + vec1.y*vec2.y) /
  262.                         SQ(distance(0.0,0.0,vec1.x,vec1.y));
  263.  
  264.         vec1.x  = l2.x2 - l2.x1;
  265.         vec1.y = l2.y2 - l2.y1;
  266.  
  267.         vec2.x = vert.d1 - l2.x1;
  268.         vec2.y = vert.d2 - l2.y1;
  269.  
  270.         if (vec2.x == 0.0 && vec2.y == 0.0) /* Zero length !! */
  271.                 parms[1] = -1.0;
  272.         else
  273.                 parms[1] = (vec1.x*vec2.x + vec1.y*vec2.y) /
  274.                         SQ(distance(0.0,0.0,vec1.x,vec1.y));
  275.  
  276. }
  277.