Graphic algorithm Unions, intersect, subtract

I need a good source for reading up on how to create a algorithm to take two polylines (a path comprised of many lines) and performing a union, subtraction, or intersection between them. This is tied to a custom API so I need to understand the underlying algorithm.

Plus any sources in a VB dialect would be doubly helpful.

From stackoverflow

• Several routines for you here. Hope you find them useful :-)

  // routine to calculate the square of either the shortest distance or largest distance
  // from the CPoint to the intersection point of a ray fired at an angle flAngle 
  // radians at an array of line segments
  // this routine returns TRUE if an intersection has been found in which case flD
  // is valid and holds the square of the distance.
  // and returns FALSE if no valid intersection was found
  // If an intersection was found, then intersectionPoint is set to the point found
  bool CalcIntersection(const CPoint &cPoint, 
                        const float flAngle,
            const int nVertexTotal,
            const CPoint *pVertexList,
            const BOOL bMin,
            float &flD,
            CPoint &intersectionPoint)
  
  {
      float d, dsx, dsy, dx, dy, lambda, mu, px, py;
      int p0x, p0y, p1x, p1y;
  
      // get source position 
      const float flSx = (float)cPoint.x;
      const float flSy = -(float)cPoint.y;
  
      // calc trig functions
      const float flTan = tanf(flAngle);
      const float flSin = sinf(flAngle);
      const float flCos = cosf(flAngle);
      const bool bUseSin = fabsf(flSin) > fabsf(flCos);
  
      // initialise distance
      flD = (bMin ? FLT_MAX : 0.0f);
  
      // for each line segment in protective feature
      for(int i = 0; i < nVertexTotal; i++) 
      {
       // get coordinates of line (negate the y value so the y-axis is upwards)
       p0x = pVertexList[i].x;
       p0y = -pVertexList[i].y;
       p1x = pVertexList[i + 1].x;
       p1y = -pVertexList[i + 1].y;
  
       // calc. deltas
       dsx = (float)(cPoint.x - p0x);
       dsy = (float)(-cPoint.y - p0y);
       dx = (float)(p1x - p0x);
       dy = (float)(p1y - p0y);
  
       // calc. denominator
       d = dy * flTan - dx;
  
       // if line & ray are parallel
       if(fabsf(d) < 1.0e-7f)
        continue;
  
       // calc. intersection point parameter
       lambda = (dsy * flTan - dsx) / d;
  
       // if intersection is not valid
       if((lambda <= 0.0f) || (lambda > 1.0f))
        continue;
  
       // if sine is bigger than cosine
       if(bUseSin){
        mu = ((float)p0x + lambda * dx - flSx) / flSin;
       } else {
        mu = ((float)p0y + lambda * dy - flSy) / flCos;
       }
  
       // if intersection is valid
       if(mu >= 0.0f){
  
        // calc. intersection point
        px = (float)p0x + lambda * dx;
        py = (float)p0y + lambda * dy;
  
        // calc. distance between intersection point & source point
        dx = px - flSx;
        dy = py - flSy;
        d = dx * dx + dy * dy;
  
        // compare with relevant value
        if(bMin){
         if(d < flD)
         {
          flD = d;
          intersectionPoint.x = RoundValue(px);
          intersectionPoint.y = -RoundValue(py);
         }
        } else {
         if(d > flD)
         {
          flD = d;
          intersectionPoint.x = RoundValue(px);
          intersectionPoint.y = -RoundValue(py);
         }
        }
       }
      }
  
      // return 
      return(bMin ? (flD != FLT_MAX) : (flD != 0.0f));
  }
  
  // Routine to calculate the square of the distance from the CPoint to the
  // intersection point of a ray fired at an angle flAngle radians at a line.
  // This routine returns TRUE if an intersection has been found in which case flD
  // is valid and holds the square of the distance.
  // Returns FALSE if no valid intersection was found.
  // If an intersection was found, then intersectionPoint is set to the point found.
  bool CalcIntersection(const CPoint &cPoint, 
                        const float flAngle,
            const CPoint &PointA,
            const CPoint &PointB,
            const bool bExtendLine,
            float &flD,
            CPoint &intersectionPoint)
  {
      // get source position 
      const float flSx = (float)cPoint.x;
      const float flSy = -(float)cPoint.y;
  
      // calc trig functions
      float flTan = tanf(flAngle);
      float flSin = sinf(flAngle);
      float flCos = cosf(flAngle);
      const bool bUseSin = fabsf(flSin) > fabsf(flCos);
  
      // get coordinates of line (negate the y value so the y-axis is upwards)
      const int p0x = PointA.x;
      const int p0y = -PointA.y;
      const int p1x = PointB.x;
      const int p1y = -PointB.y;
  
      // calc. deltas
      const float dsx = (float)(cPoint.x - p0x);
      const float dsy = (float)(-cPoint.y - p0y);
      float dx = (float)(p1x - p0x);
      float dy = (float)(p1y - p0y);
  
      // Calc. denominator
      const float d = dy * flTan - dx;
  
      // If line & ray are parallel
      if(fabsf(d) < 1.0e-7f)
       return false;
  
      // calc. intersection point parameter
      const float lambda = (dsy * flTan - dsx) / d;
  
      // If extending line to meet point, don't check for ray missing line
      if(!bExtendLine)
      {
       // If intersection is not valid
       if((lambda <= 0.0f) || (lambda > 1.0f))
        return false; // Ray missed line
      }
  
      // If sine is bigger than cosine
      float mu;
      if(bUseSin){
       mu = ((float)p0x + lambda * dx - flSx) / flSin;
      } else {
       mu = ((float)p0y + lambda * dy - flSy) / flCos;
      }
  
      // if intersection is valid
      if(mu >= 0.0f)
      {
       // calc. intersection point
       const float px = (float)p0x + lambda * dx;
       const float py = (float)p0y + lambda * dy;
  
       // calc. distance between intersection point & source point
       dx = px - flSx;
       dy = py - flSy;
       flD = (dx * dx) + (dy * dy);
  
       intersectionPoint.x = RoundValue(px);
       intersectionPoint.y = -RoundValue(py);
       return true;
      }
  
      return false;
  }
  
  // Fillet (with a radius of 0) two lines. From point source fired at angle (radians) to line Line1A, Line1B.
  // Modifies line end point Line1B. If the ray does not intersect line, then it is rotates every 90 degrees
  // and tried again until fillet is complete.
  void Fillet(const CPoint &source, const float fThetaRadians, const CPoint &Line1A, CPoint &Line1B)
  {
      if(Line1A == Line1B)
       return; // No line
  
      float dist;
  
      if(CalcIntersection(source, fThetaRadians, Line1A, Line1B, true, dist, Line1B))
       return;
      if(CalcIntersection(source, CalcBaseFloat(TWO_PI, fThetaRadians + PI * 0.5f), Line1A, Line1B, true, dist, Line1B))
       return;
      if(CalcIntersection(source, CalcBaseFloat(TWO_PI, fThetaRadians + PI), Line1A, Line1B, true, dist, Line1B))
       return;
      if(!CalcIntersection(source, CalcBaseFloat(TWO_PI, fThetaRadians + PI * 1.5f), Line1A, Line1B, true, dist, Line1B))
       ASSERT(FALSE); // Could not find intersection?
  }
  
  // routine to determine if an array of line segments cross gridSquare
  // x and y give the float coordinates of the corners
  BOOL CrossGridSquare(int nV, const CPoint *pV, 
                       const CRect &extent, const CRect  &gridSquare)
  {
      // test extents
      if( (extent.right < gridSquare.left) ||
       (extent.left > gridSquare.right) ||  
       (extent.top  > gridSquare.bottom) ||
       (extent.bottom < gridSquare.top))
      {
       return FALSE;
      }
  
      float a, b, c, dx, dy, s, x[4], y[4];
      int max_x, max_y, min_x, min_y, p0x, p0y, p1x, p1y, sign, sign_old;
  
      // construct array of vertices for grid square
      x[0] = (float)gridSquare.left;
      y[0] = (float)gridSquare.top;
      x[1] = (float)(gridSquare.right);
      y[1] = y[0];
      x[2] = x[1];
      y[2] = (float)(gridSquare.bottom);
      x[3] = x[0];
      y[3] = y[2];
  
      // for each line segment
      for(int i = 0; i < nV; i++) 
      {
       // get end-points
       p0x = pV[i].x;
       p0y = pV[i].y;
       p1x = pV[i + 1].x;
       p1y = pV[i + 1].y;
  
       // determine line extent
       if(p0x > p1x){
        min_x = p1x;
        max_x = p0x;
       } else {
        min_x = p0x;
        max_x = p1x;
       }
  
       if(p0y > p1y){
        min_y = p1y;
        max_y = p0y;
       } else {
        min_y = p0y;
        max_y = p1y;
       }
  
       // test to see if grid square is outside of line segment extent
       if( (max_x < gridSquare.left)  ||
        (min_x > gridSquare.right) ||  
        (max_y < gridSquare.top)   ||
        (min_y > gridSquare.bottom))
       {
        continue;
       }
  
       // calc. line equation
       dx = (float)(p1x - p0x);
       dy = (float)(p1y - p0y);
       a = dy;
       b = -dx;
       c = -dy * (float)p0x + dx * (float)p0y;
  
       // evaluate line eqn. at first grid square vertex
       s = a * x[0] + b * y[0] + c;
       if(s < 0.0f){
        sign_old = -1;
       } else if(s > 1.0f){
        sign_old = 1;
       } else {
        sign_old = 0;
       }
  
       // evaluate line eqn. at other grid square vertices
       for (int j = 1; j < 4; j++) 
       {
        s = a * x[j] + b * y[j] + c;
        if(s < 0.0f){
         sign = -1;
        } else if(s > 1.0f){
         sign = 1;
        } else {
         sign = 0;
        }
  
        // if there has been a chnage in sign
        if(sign != sign_old)
         return TRUE;
       }
      }
  
      return FALSE;
  }
  
  // calculate the square of the shortest distance from point s
  // and the line segment between p0 and p1
  // t is the point on the line from which the minimum distance
  // is measured
  float CalcShortestDistanceSqr(const CPoint &s,
                 const CPoint &p0,
                 const CPoint &p1,
              CPoint &t)
  {
      // if point is at a vertex
      if((s == p0) || (s == p1))
       return(0.0F);
  
      // calc. deltas
      int dx = p1.x - p0.x;
      int dy = p1.y - p0.y;
      int dsx = s.x - p0.x;
      int dsy = s.y - p0.y;
  
      // if both deltas are zero 
      if((dx == 0) && (dy == 0))
      {
       // shortest distance is distance is to either vertex
       float l = (float)(dsx * dsx + dsy * dsy);
       t = p0;
       return(l);
      }
  
      // calc. point, p, on line that is closest to sourcePosition
      // p = p0 + l * (p1 - p0)
      float l = (float)(dsx * dx + dsy * dy) / (float)(dx * dx + dy * dy);
  
      // if intersection is beyond p0
      if(l <= 0.0F){
  
       // shortest distance is to p0
       l = (float)(dsx * dsx + dsy * dsy);
       t = p0;
  
      // else if intersection is beyond p1
      } else if(l >= 1.0F){
  
       // shortest distance is to p1
       dsx = s.x - p1.x;
       dsy = s.y - p1.y;
       l = (float)(dsx * dsx + dsy * dsy);
       t = p1;
  
      // if intersection is between line end points
      } else {
       // calc. perpendicular distance
       float ldx = (float)dsx - l * (float)dx;
       float ldy = (float)dsy - l * (float)dy;
       t.x = p0.x + RoundValue(l * (float)dx);
       t.y = p0.y + RoundValue(l * (float)dy);
       l = ldx * ldx + ldy * ldy;
      }
  
      return(l);
  }
  
  // Calculates the bounding rectangle around a set of points
  // Returns TRUE if the rectangle is not empty (has area), FALSE otherwise
  // Opposite of CreateRectPoints()
  BOOL CalcBoundingRectangle(const CPoint *pVertexList, const int nVertexTotal, CRect &rect)
  {
      rect.SetRectEmpty();
      if(nVertexTotal < 2)
      {
       ASSERT(FALSE); // Must have at least 2 points
       return FALSE;
      }
  
      // First point, set rectangle (no area at this point)
      rect.left = rect.right = pVertexList[0].x;
      rect.top = rect.bottom = pVertexList[0].y;
  
      // Increst rectangle by looking at other points
      for(int n = 1; n < nVertexTotal; n++)
      {
       if(rect.left > pVertexList[n].x) // Take minimum
        rect.left = pVertexList[n].x;
  
       if(rect.right < pVertexList[n].x) // Take maximum
        rect.right = pVertexList[n].x;
  
       if(rect.top > pVertexList[n].y)  // Take minimum
        rect.top = pVertexList[n].y;
  
       if(rect.bottom < pVertexList[n].y) // Take maximum
        rect.bottom = pVertexList[n].y;
      }
  
      rect.NormalizeRect(); // Normalise rectangle
      return !(rect.IsRectEmpty());
  }
  

• This catalogue of implementations of intersection algorithms from the Stony Brook Algorithm Repository might be useful. The repository is managed by Steven Skiena, author of a very well respected book on algorithms: The Algorithm Design Manual.

  That's his own Amazon exec link by the way :)

  RS Conley : That looks like it will do the trick and just a important I know what the algorithm is called and found more stuff on google, thanks.

  MarkJ : Knowing the name of the algorithm is half the battle - maybe more.