// Example solution to programming assignment 7.
// CS 509, Spring 2002.

#include <algorithm>
#include <cassert>
#include <cstdio>
#include <functional>
#include <iostream>
#include <map>
#include <string>
#include <sys/time.h>
#include <unistd.h>
#include <vector>
#include <strings.h>

typedef std::pair<int, int> point;
typedef std::vector<point>  point_vec;

//
// This is the slow, O(n^3) solution for finding the largest collinear point
// set.
//
// The algorithm is straightforward.  Each pair of points in the set form a
// line.  Checking every other point in the set against the given line
// determines how many points are collinear to that line (there must be at
// least two such points.  Doing that for every pair of points in the set and
// remembering the line containing the largest number of points determines a
// maximum collinear point set.
//
// There O(n) points to compare to a line, and there are O(n^2) pairs of points
// in the set.  Doing O(n) work O(n^2) times gives an O(n^3) algorithm
//
// The innermost loop should be replace by a std::for_each() call with a
// functior.
//

static point_vec 
scheck_colinear_set(const point_vec & points) {

  point_vec pts;

  if (points.size() < 2)
    return pts;

  // Invariant: pts contains a maximum set of collinear points that contain any
  // of points[0..i-1].
  
  for (unsigned i = 0; i < points.size(); i++) {

    // Invariant: new-set contains a maxium set of collinear points that
    // contain any of the points points[i..j-1].

    for (unsigned j = i + 1; j < points.size(); j++) {

      int x_step = points[j].first - points[i].first;
      int y_step = points[j].second - points[i].second;

      if ((x_step == 0) and (y_step == 0)) {
	std::cerr << "!!! The input point set has duplicates of the point "
	          << points[j] << ".\n";
	exit(EXIT_FAILURE);
	}

      if (x_step == 0) 
	y_step = 1;
      else if (y_step == 0)
	x_step = 1;
      else {
	const int c = gcd(x_step, y_step);
	x_step /= c;
	y_step /= c;
	}

      point_vec new_set;
      new_set.push_back(points[i]);
      new_set.push_back(points[j]);

      // Invariant: new-set contains a maxium set of collinear points that
      // contain any of the points points[i..k-1].

      for (unsigned k = j + 1; k < points.size(); k++) {
	int delta_x =  points[k].first - points[i].first;
	int delta_y =  points[k].second - points[i].second;
	if (delta_x == 0)
	  delta_y = 1;
	if (delta_y == 0)
	  delta_x = 1;
	if (same_step(x_step, y_step, delta_x, delta_y))
	  new_set.push_back(points[k]);
	}

      if (new_set.size() > pts.size())
	pts = new_set;
      }
    }

  return pts;
  }

//
// We need to do a little housekeeping before we can get to a faster algorithm
// for finding maximum collinear-point sets.  
// 
// Representing slopes is a problem that requires some attention.  Floating
// point numbers seem like a good choice, but they have two difficulties.
// First, some lines have inifinite slopes, which is difficult to handle with
// floating-point numbers (actually not that difficult with the IEEE standard).
// The second problem with floating-point numbers is accuracy; it's possible
// that two points that are collinear with a third may have two different
// slopes if they're computed with floating-point numbers, particularly if the
// line is almost vertical or horizontal.
//
// The easiest way to represent slope with accuracy is to do everything but
// the final divide to get a single number: a slope is represented as a
// rational number with the numerator representing the delta x and the
// denominator representing the delta y.
//

class rational {

  public:

    rational(int num, int den) { 

      if (den == 0)
	s = infinite;
      else if (num == 0)
	s = zero;
      else if (((num < 0) and (den < 0)) or ((num > 0) and (den > 0)))
	s = positive;
      else
	s = negative;

      n = num < 0 ? -num : num;
      d = den < 0 ? -den : den;

      if ((d == 0) and (n == 0)) {
	std::cerr << "!!! Duplicate points found in input.\n";
	exit(EXIT_FAILURE);
	}
      else if (d == 0) 
	n = 1;
      else if (n == 0)
	d = 1;
      else {
	const int c = gcd(d, n);
	d /= c;
	n /= c;
	}
      }

    bool operator < (const rational & r) const { 
      if (r.s == infinite) 
	return s != infinite;
      else
	return (s < r.s) or 
	       ((s == r.s) and (n*r.d < r.n*d)); 
      }

  private:

    enum sign {
      negative,
      zero,
      positive,
      infinite
      };

    sign s;
    int n;
    int d;
  };

//
// A bit more housekeeping.  A map indexed by rationals should give you a clue
// about how the faster algorithm's going to work.
//

typedef std::map<rational, point_vec> slope_list;
typedef slope_list::iterator          slist_iter;
typedef slope_list::const_iterator    slist_citer;
typedef slope_list::value_type        slist_etype;

//
// Just a little more housekeeping.  Return true if and only if the point
// vector in the first slist map element is smaller than the point vector in
// the second slist map element.
//

static bool
operator < (const slist_etype & p1, const slist_etype & p2) {
  return p1.second.size() < p2.second.size();
  }

//
// Ok, now.  This is the faster O(n^2 log n) solution for finding the largest
// collinear point set.  The trick is to pick a point and find the slopes
// between that point and all the other points in the set; points collinear
// with the chosen point will all have the same slope.  Storing the points in a
// map indexed by slope is a simple and relatively efficient way of collecting
// points of identical slope.  Each map access takes O(log n) and the map will
// be accessed O(n) times for each chosen point.  There are O(n) chosen points,
// so the algorithm has O(n)O(n)O(log n) = O(n^2 log n) running time, as
// advertised.
//

static point_vec
fcheck_colinear_set(const point_vec & points) {

  point_vec pts;

  if (points.size() < 2)
    return pts;

  // Invariant: pts contains a maximum collinear set of points that contain
  // points from points[0..i-1].

  for (unsigned i = 0; i < points.size() - 1; i++) {

    slope_list slopes;
    const point & pt = points[i];

    // Invariant: the points points[i+1..j-] have been placed in a line that
    // contains the point points[i].

    for (unsigned j = i + 1; j < points.size(); j++)
      slopes[rational(points[j].first - pt.first, 
		      points[j].second - pt.second)].push_back(points[j]);

    const slist_iter e = slopes.end();
    const slist_citer i = std::max_element(slopes.begin(), e);
    assert(i != e);
    
    if (i->second.size() >= pts.size()) {
      pts = i->second;
      pts.push_back(pt);
      }
    }

  return pts;
  }

//
// Deja vu, c++ style.
//

extern point_vec colinear_set(const point_vec &);
extern int gcd(int, int);
extern std::ostream & operator << (std::ostream &, const point &);
extern std::istream & operator >> (std::istream &, const point &);

static void check_times(time_t, time_t, unsigned, unsigned, const char * const);
static void run_test(bool);
static point_vec read_points(std::istream &);
static void check_colinear(const point_vec &, const std::string &);
static point_vec scheck_colinear_set(const point_vec &);
static point_vec fcheck_colinear_set(const point_vec &);
static void print_time(time_t);
static point_vec timed_run(point_vec (*)(const point_vec &), const point_vec &, time_t &);

//
// A functor that determines if the point passed in via operator() is colinear
// with the two points passed in when the functor was constructed.
//

struct
noncolinear_helper {

  // Because points have integer coordinates, a point is collinear with a line
  // if it is an integer number of steps away from some designated origin point
  // on a line.

  const point origin;
  int x_step, y_step;


  // Construct a noncolinear helper to determine if points are collinear with
  // the line defined by p1 and p2.

  noncolinear_helper(const point & p1, const point & p2) : origin(p1){

    x_step = origin.first - p2.first;
    y_step = origin.second - p2.second;

    assert(x_step or y_step);
    
    if (x_step == 0)
      y_step = 1;
    else if (y_step == 0)
      x_step = 1;
    else {
      const int c = gcd(x_step, y_step);
      x_step /= c;
      y_step /= c;
      }
    }


  // Return true iff p is not collinear with the line defined when this
  // helper was constructed.

  bool operator()(point p) {

    const int 
      delta_x = origin.first - p.first,
      delta_y = origin.second - p.second,
      i = x_step ? delta_x/x_step : delta_y/y_step;

    return ((i*x_step != delta_x) or (i*y_step != delta_y));
    }
  };

//
// Determine if the points in pts are collinear.
//

static void
check_colinear(const point_vec & pts, const std::string & what) {

  if (pts.size() > 2) {
    point_vec::const_iterator e = pts.end();
    point_vec::const_iterator f =
      std::find_if(pts.begin() + 2, e, noncolinear_helper(pts[0], pts[1]));

    if (f != e) {
      std::cerr << "!!! Non-colinear points " << pts[0] << ", " << pts[1]
	        << ", " << *f << "\n!!! found in " << what << ".\n";
      exit(EXIT_FAILURE);
      }
    }
  }

//
// Compare an answer time to a check time and print a message if the answer
// time's too slow compared to the check time.
//

static void
check_times(
  time_t answert, time_t checkt,
  unsigned asize, unsigned ssize, 
  const char * const what) { 

  if (answert > 2*checkt) {
    fprintf(stderr, "!!! %s ", what);
    print_time(answert);
    fprintf(stderr, " vs. ");
    print_time(checkt);
    fprintf(stderr, " (x%3.1f) ",
	    static_cast<double>(answert)/static_cast<double>(checkt));
    fprintf(stderr, "for %d of %d points.\n", asize, ssize);
    exit(EXIT_FAILURE);
    }
  }

//
// A function that needs no introduction.
//

int
gcd(int a, int b) {

  assert((a > 0) and (b > 0));

  while (a != b)
    if (a > b) {
      int t = a;
      a = b;
      b = t - b;
      }
    else {
      int t = b;
      b = a;
      a = t - a;
      }

  return a;
  }

//
// The command-line argument determines whether or not the times to find the
// maximum collinear-point sets are printed or not.
//

int
main(int argc, char * argv[]) {
  
  if (argc > 2) {
    std::cerr << "Command format is \"" << argv[0] << " [-t]\".\n";
    return EXIT_FAILURE;
    }

  bool telltime = false;

  if (argc == 2)
    if (strcmp(argv[1], "-t") == 0)
      telltime = true;
    else {
      std::cerr << "Command format is \"" << argv[0] << " [-t]\".\n";
      return EXIT_FAILURE;
      }

  run_test(telltime);
  }

//
// Print a point to an output stream.
//

std::ostream & 
operator << (std::ostream & os, const point & p) {
  return os << "(" << p.first << ", " << p.second << ")";
  }

//
// Read a point from an input stream; die with a message if there's an error.
//

std::istream & 
operator >> (std::istream & is, point & p) {

  is >> p.first;
  if (is.eof()) return is;
  if (!is.good()) {
    std::cerr << "Input error during x-coordinate input.\n";
    exit(EXIT_FAILURE);
    }

  if (!(is >> p.second)) {
    std::cerr << "Error during y-coordinate input.\n";
    exit(EXIT_FAILURE);
    }

  return is;
  }

//
// Print a time in microseconds using a small number of digits by changing the
// time scale.
//

static void 
print_time(time_t t) {
  if (t < 1000)
    fprintf(stderr, "%5u usec", static_cast<unsigned>(t));
  else if (t < 100000)
    fprintf(stderr, "%5.1f msec", static_cast<double>(t)/1000.0);
  else
    fprintf(stderr, "%5.1f  sec", static_cast<double>(t)/1000000.0);
  }

//
// Read the points found in the input stream and return them in a point vector.
// This seems like an ideal place to use an input-stream iterator, but g++
// woun't cooperate by compiling it.l
//

static point_vec
read_points(std::istream & ins) {

  point_vec pts;
  point pt;

  while (std::cin >> pt) 
    pts.push_back(pt);

  return pts;
  }

//
// Run the external colinear() procedure against the two defined here and
// compare the results for correctness and speed.
//

static void
run_test(bool telltimes) {

  time_t answert, fast_et, slow_et;
  const point_vec pts = read_points(std::cin);

  const point_vec ans = timed_run(colinear_set, pts, answert);
  check_colinear(ans, "your colinear_set() answer");

  const point_vec my_sans = timed_run(scheck_colinear_set, pts, slow_et);
  check_colinear(my_sans, "the slow-check answer");

  const point_vec my_fans = timed_run(fcheck_colinear_set, pts, fast_et);
  check_colinear(my_fans, "the fast-check answer");

  if (ans.size() < my_sans.size()) {
    std::cerr << "!!! The computed colinear-point set size " << ans.size() 
              << " is smaller than the check colinear-point set size " 
	      << my_sans.size() << ".\n";
    exit(EXIT_FAILURE);
    }

  if (my_sans.size() < ans.size()) {
    std::cerr << "!!! The check colinear-point set size " << my_sans.size() 
              << " is smaller than the computed colinear-point set size " 
	      << ans.size() << ".\n";
    exit(EXIT_FAILURE);
    }

  if (my_sans.size() != my_fans.size()) {
    std::cerr << "!!! The slow-check colinear-point set size " << my_sans.size() 
              << " is not equal to the fast check colinear-point set size " 
	      << my_fans.size() << ".\n";
    exit(EXIT_FAILURE);
    }

  if (telltimes) {
    fprintf(stderr, "Your colinear_set() took ");
    print_time(answert);
    fprintf(stderr, " to find an answer of\n");
    fprintf(stderr, "%d points from a set of %d points.\n", 
	    ans.size(), pts.size());
    fprintf(stderr, "The slow colinear_set() check took ");
    print_time(slow_et);
    fprintf(stderr, " to find an answer.\n");
    fprintf(stderr, "The fast colinear_set() check took ");
    print_time(fast_et);
    fprintf(stderr, " to find an answer.\n");
    }

  check_times(answert, slow_et, ans.size(), pts.size(), "slow");
  check_times(answert, fast_et, ans.size(), pts.size(), "fast");
  }

//
// Run the collinear procedure cs against the clock using the point set pts and
// returning the maximum collinear-point set found.  et contains the ellapsed
// time in microseconds it took to find the set.
//

static point_vec
timed_run(
  point_vec (* cs)(const point_vec &), const point_vec & pts, time_t & et) {

  struct timeval s, e;

  gettimeofday(&s, NULL);
  point_vec cps = cs(pts);
  gettimeofday(&e, NULL);

  if (s.tv_sec < e.tv_sec) { 
    e.tv_sec--;
    e.tv_usec += 1000000; 
    }
  et = (e.tv_sec - s.tv_sec)*1000000 + (e.tv_usec - s.tv_usec);

  return cps;
  }

// 
// And that's it.
//


syntax highlighted by Code2HTML, v. 0.9