Algorithm Alley
by Thomas Gettys

Listing One
//***********************************************************************
// This is the actual Graham scan algorithm. It is called after the pivot 
// point has been located and the rest of the points have been sorted by 
// their polar angle.

void Graham_Scan(point S[], int N)
{
   int i;
   point p;

   Push(S[0]);
   Push(S[1]);
   Push(S[2]);

   for (i = 3; i < N; i++)
   {
      while (CrossProduct(Stack2nd, StackTop, S[i]) < 0)
      {
         Pop(&p);
      };
      Push(S[i]);
   }
}


Listing Two
typedef struct
{
   int x, y;
} point;

// Find point with smallest y coordinate; if there is more
// than one select the one with the smallest x coordinate.

int LocatePivot(point S[], int N)
{
   int i, ndx, minx, miny;

// initialize minimum coordinates found so far
   ndx  = 0;
   minx = S[ndx].x;
   miny = S[ndx].y;

// locate lowest point in set
   for (i = 1; i < N; i++)
   {
      if (S[i].y > miny) continue;

//    found a candidate for a new pivot
      if (S[i].y < miny)
      {
         ndx = i;
         miny = S[i].y;
      }
      else if (S[i].x < minx)
      {
         ndx = i;
         minx = S[i].x;
      }
   }

// all done - return index of pivot point
   return(ndx);
}

Listing Three
#define sizeof_Stack 100
#define StackTop     stack[TOS_ptr-1]
#define Stack2nd     stack[TOS_ptr-2]

int TOS_ptr = 0;           // Top Of Stack pointer
point stack[sizeof_Stack]; // the actual point stack array

//*********************************************************
// pushes point p onto the stack
// if successful Push() returns 1
// on stack overflow Push() returns 0

int Push(point p)
{
   if (TOS_ptr == sizeof_Stack) return(0);

   stack[TOS_ptr++] = p;
   return(1);
}
//*********************************************************
// pops a point off of the stack
// if successful Pop() returns 1
// on stack underflow Pop() returns 0

int Pop(point *p)
{
   if (TOS_ptr == 0) return(0);
   *p = stack[--TOS_ptr];
   return(1);
}

Listing Four
//*************************************************************************
// Computes the cross product (p1-p0)x(p2-p0). If p2 is to the left of the 
// line from p0 to p1 the result will be positive.

int CrossProduct(point p0, point p1, point p2)
{
   int xprod;
   xprod = (p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y);
   return(xprod);
}
//*********************************************************
// Support function for qsort()
// a < b if a makes a smaller angle with Pivot than does b

int compare_points(const void *a, const void *b)
{
   int xprod;
   point *point_a, *point_b;

   point_a = (point *)a;
   point_b = (point *)b;

   xprod = CrossProduct(Pivot, *point_a, *point_b);

   if (xprod < 0) return(+1);    // a > b
   else
   if (xprod > 0) return(-1);    // a < b

   return(0);                    // a = b
}

Listing Five
point Pivot, S[100];

void main(void)
{
   int i, N;
   point p;

// make up some points
   N = 12;
   randomize();
   for (i = 0; i < N; i++)
   {
      S[i].x = random(10);
      S[i].y = random(10);
      printf("(%3d, %3d)\n", S[i].x, S[i].y);
   }
// find the convex hull
   i = LocatePivot(S, N);
   Pivot = S[i];
   S[i]  = S[0];
   S[0]  = Pivot;

   qsort((void *)&S[1], N-1, sizeof(S[0]), compare_points);

   Graham_Scan(S, N);

// output points of convex hull in clockwise order
   for (i = 0; Pop(&p) != 0; i++)
   {
      printf("(%3d, %3d)\n", p.x, p.y);
   }
   printf("%d points in convex hull\n", i);
}




3