Friday, December 7, 2012

MCPC 2012 - Problem D : Minimum bounding rectangle

May Allah`s peace , mercy and blessing be upon you

P.S: This is not the original text, I rephrased and summarized the original one because I'm lazy and I couldn't type the whole text. if there is an error it's my mistake.

[D] SCPC:
We need to build a contest hall in the shape of a cuboid (3D rectangle). for economical reasons, we need to minimize the size (area) of the rectangle. So given the coordinates of the teams, calculate the minimum rectangle area of land for a building big enough to fit all the teams.

Input :
T: number of cases.
N: number of teams.
Xi Yi: coordinates of each team.
1<=T,N<=100,
-10000<=Xi,Yi<=10000

Output:
Single line starting with case number, followed by the area of the minimum rectangle rounded to 4 decimal places.

Example :
4
4
0 0
1 1
1 0
0 1
4
0 0
2 2
2 0
0 2
4
0 0
1 2
1 0
0 2
3
0 0
1 0
1 1
--->
Case 1: 1.0000
Case 2: 4.0000
Case 3: 2.0000
Case 4: 1.0000

Contest Analysis :
Actually I didn't read this problem during the competition, so this solution wasn't tested against judge's inputs and therefore it's not conclusive.
This is a classic problem called minimum bounding rectangle, minimum bounding envelope and in 3D called minimum bounding box. To achieve this I am going to use the convex hull, because there is a theorem that says "The minimum bounding rectangle's orientations are determined by the polygon edges". So the steps are :
  • Compute the convex hull of the cloud.
  • For each edge of the convex hull :
    • Compute the edge orientation.
    • Rotate the convex hull using this orientation in order to compute easily the bounding rectangle area.
    • Store the orientation corresponding to the minimum area found.
  • Return the rectangle corresponding to the minimum area found.

Solution in C++:
A very long code in C++ , during a contest this should be done in Java.
To find the convex hull, I translated a Javascript code described by Teemu Korhonen in this blog post.

#include <iostream>
#include <iomanip>
#include <fstream>
#include <vector>
#include <algorithm>
#define  _USE_MATH_DEFINES
#include <cmath>
#include <stack>

#define For(i,n) for(int i(0),_n(n);i<_n;++i)

using namespace std;

class Point {
  public:
 double x;
 double y;
 double angle;
};

bool comp(double a,double b) {return (a<b);}
bool compy(Point p1,Point p2) {return (p1.y<p2.y);}
bool compx(Point p1,Point p2) {return (p1.x<p2.x);}
bool compangle(Point p1,Point p2) {return (p1.angle<p2.angle);}

void rotate(Point &point, Point pivot, double angle) {
 float s = sin(angle);
 float c = cos(angle);
 point.x -= pivot.x;
 point.y -= pivot.y;
 double newx = point.x*c - point.y*s;
 double newy = point.x*s + point.y*c;
 point.x = newx + pivot.x;
 point.y = newy + pivot.y;
}

double ptdistance(Point p1,Point p2) {
 return pow(p2.x-p1.x,2) + pow(p2.y-p1.y,2);
}

double calculate_angle(Point point, Point pivot) {
 double x = point.x - pivot.x;
 double y = point.y - pivot.y;
 if(x>0 && y>0) return atan(y/x);
 if(x<0 && y>0) return atan(-x/y) + M_PI/2;
 if(x<0 && y<0) return atan(y/x) + M_PI;
 if(x>0 && y<0) return atan(-x/y) + M_PI/2 + M_PI;
 if(x==0 && y>0) return M_PI/2;
 if(x<0 && y==0) return M_PI;
 if(x==0 && y<0) return M_PI/2 + M_PI;
 else return 0;
}

// Three points are a counter-clockwise turn if ccw > 0, clockwise if ccw < 0, and collinear if ccw = 0 
double ccw(Point p1, Point p2, Point p3) { 
 return (p2.x - p1.x)*(p3.y - p1.y) - (p2.y - p1.y)*(p3.x - p1.x); 
}

double parallelsurface(vector<Point> hull, Point point, Point pivot) {
 double angle = calculate_angle(point,pivot);
 vector<Point> newhull(hull);
 for(vector<Point>::iterator it=newhull.begin(); it!=newhull.end(); ++it) {
  if(angle < M_PI/2) rotate(*it,pivot,-angle);
  if(angle > M_PI/2 && angle< M_PI) rotate(*it,pivot,-angle+M_PI/2);
  if(angle > M_PI && angle< 3*M_PI/2) rotate(*it,pivot,-angle+M_PI);
  if(angle > 3*M_PI/2 && angle< 2*M_PI) rotate(*it,pivot,-angle+3*M_PI/2);
 }
 double minx = (*min_element(newhull.begin(), newhull.end(), compx)).x;
 double miny = (*min_element(newhull.begin(), newhull.end(), compy)).y;
 double maxx = (*max_element(newhull.begin(), newhull.end(), compx)).x;
 double maxy = (*max_element(newhull.begin(), newhull.end(), compy)).y;
 return (maxx-minx)*(maxy-miny);
}

int main(int argc, char const *argv[]) {
 long long T;
 fstream f("scpc.in");
 f >> T;
 For(t,T) {
  long long N;
  vector<Point> points;
  vector<Point> hull;
  vector<double> surfaces;
  f >> N;
  For(i,N) {
   Point pt;
   f >> pt.x;
   f >> pt.y;
   points.push_back(pt);
  }
  // All points with positive coordinates
  int minx = (*min_element(points.begin(), points.end(), compx)).x;
  int miny = (*min_element(points.begin(), points.end(), compy)).y;
  for(vector<Point>::iterator it=points.begin(); it!=points.end(); ++it) {
   if(minx<0) (*it).x += abs(minx);
   if(miny<0) (*it).y += abs(miny);
  }

  // Get reference point C (lowest leftmost) and remove it from points
  vector<Point>::iterator ref = points.begin();
  for(vector<Point>::iterator it=points.begin(); it!=points.end(); ++it) {
   if((*it).y < (*ref).y) ref = it;
   if((*it).y == (*ref).y)
    if((*it).x > (*ref).x) ref = it;
  }
  Point C;
  C.x = (*ref).x;
  C.y = (*ref).y;
  C.angle = 0;
  points.erase(ref);
  --N;

  // Stack to use for removing points to not interfere with the loop
  stack<vector<Point>::iterator> rmpts; 

  // Fill the angle field and 
  for(vector<Point>::iterator it=points.begin(); it!=points.end(); ++it) {
   (*it).angle = calculate_angle(*it,C);
   if((*it).x == C.x && (*it).y == C.y) {
    rmpts.push(it);
    --N;
   }
  }
  while(!rmpts.empty()) { // Effective removal
   points.erase(rmpts.top());
   rmpts.pop();
  }

  // Sort by angle
  sort(points.begin(), points.end(), compangle);

  // If two points have same angle remove the closest to C
  for(vector<Point>::iterator it=points.begin()+1, prev=points.begin(); it!=points.end(); prev=it, ++it) {
   if((*prev).angle == (*it).angle) {
    double d1 = ptdistance(C, *it);
    double d2 = ptdistance(C, *prev);
    if(d1 > d2) rmpts.push(prev);
    else rmpts.push(it);
    // I don't know why I keep updating N even if I'm done with it !!
    --N;
   }
  }
  while(!rmpts.empty()) { // Effective removal
   points.erase(rmpts.top());
   rmpts.pop();
  }

  // Building the convex hull
  hull.push_back(C);
  hull.push_back(points[0]);
  hull.push_back(points[1]);
  for(vector<Point>::iterator it=points.begin()+2; it!=points.end(); ++it) {
   Point pi = *it;
   while(ccw(hull[hull.size()-2], hull[hull.size()-1], pi) <= 0) {
    hull.pop_back();
   }
   hull.push_back(pi);
  }

  // Calculate MBR area for each side
  for(vector<Point>::iterator it=hull.begin()+1, prev=hull.begin(); it!=hull.end(); prev=it, ++it) {
   surfaces.push_back(parallelsurface(hull,*it,*prev));
  }
  surfaces.push_back(parallelsurface(hull,hull[hull.size()-1],hull[0]));

  cout << "Case " << t+1 << ": " << fixed << setprecision(4) << (*min_element(surfaces.begin(), surfaces.end(), comp)) << endl;
 }
 f.close();
 return 0;
}

No comments:

Post a Comment


Zoli Issam's blog proudly powered by Blogger | All rights reserved Jaw,B 2010-2013