ACM: 137 - Polygons

寫完132，打鐵趁熱繼續寫geometry的題目。這題很複雜，但要解決的問題很清楚，就是求兩個凸多邊形的非交集的面積(xor, 異域)。

我的算法如下：
A: 第一個凸多邊形點集合
B: 第二個凸多邊形點集合
C: A,B交集部分的點集合(初始為空集合)

1. 檢查A的每一個點，若落在B內部或邊上即將該點加入C。
2. 檢查B的每一個點，若落在A內部或邊上即將該點加入C。
3. 檢查A的每一邊，對B的每一邊逐邊檢查是否可以相交，若可相交即求交點，將交點加入C。
4. A,B交集求畢，將C集合的點逆時針排序。
5. 輸出異域面積為 (A + B – C) – C = A + B - 2C

實作方法

檢查P點是否落在圖形內部或邊上
逆時針檢查圖形的每一點，若「每一點到下一點的向量」與「每一點至P點的向量」皆為逆時針旋轉(外積>0)，則P點在圖型內部。當外積=0時，P點在圖形邊上。

檢查線段P1P2與P3P4是否相交
單純用外積檢查兩線段是否互相跨立。
當 (P1P2 x P1P3)*(P1P2 x P1P4) < 0 時代表P3, P4分別在線段P1P2的兩端，故外積異號。
當 (P3P4 x P3P1)*(P3P4 x P3P2) < 0 時代表P1, P2分別在線段P3P4的兩端，故外積異號。
當兩者皆成立時，線段P1P2與線段P3P4相交。

相交線段求交點
參考「平面內兩條線段的位置關係判定及交點求解」，用向量法求解，只需用到一次除法，減少精確度的問題。

點逆時針排序
先找最左下方的點，其他的點對該點做逆時針的極角排序。

凸多邊形面積
 
起點與終點必須為同一個點，即X0Y0 = XnYn，而所有的點必須照順時針或逆時針排序，當點為順時針排序時，總和為負數，逆時針排序時，總和為正數，但兩者絕對值相同。

TIPS
1. 製造測資時可考慮相離、相鄰、相交的情況，而相交的圖形可考慮「其中一圖有點在另一圖中」和「兩個圖形完全沒有點在另一個圖中，但有相交部份」的兩種情況。
2. 輸出時雖然每筆輸出之間不需換行，但最後結尾仍需要換行符號，否則會給WA，不是PE。

REFERENCES
平面內兩條線段的位置關係判定及交點求解

CODE

#include <iostream>
#include <vector>
using namespace std;

struct vertex{
  float x, y;
};

//return outer product of pa->pb
float cross(vertex p, vertex a, vertex b){
  return ( (a.x - p.x)*(b.y - p.y) - (b.x - p.x)*(a.y - p.y) );
}

float dist(vertex a, vertex b){
  return (b.x - a.x)*(b.x - a.x) + (b.y - a.y)*(b.y - a.y);
}

//return true if p is inside or on boundary of g
bool inside(vertex p, vector<vertex> g, int n){
  int i = n - 1;
  for(; i>=0 && (cross(g[i], g[(i+n-1)%n], p) >= 0); --i);
  return (i == -1);
}

bool intersects(vertex p1, vertex p2, vertex p3, vertex p4){
  return ((cross(p1, p2, p3) * cross(p1, p2, p4) < 0) && (cross(p3, p4, p1) * cross(p3, p4, p2) < 0));
}

vertex intersection(vertex p1, vertex p2, vertex p3, vertex p4){
  float b1 = (p2.y - p1.y)*p1.x + (p1.x - p2.x)*p1.y;
  float b2 = (p4.y - p3.y)*p3.x + (p3.x - p4.x)*p3.y;
  float d = (p2.x - p1.x)*(p4.y - p3.y) - (p4.x - p3.x)*(p2.y - p1.y);
  float d1 = b2*(p2.x - p1.x) - b1*(p4.x - p3.x);
  float d2 = b2*(p2.y - p1.y) - b1*(p4.y - p3.y);
  vertex p0;
  p0.x = d1/d;
  p0.y = d2/d;
  return p0;
}

void check_intersect(vertex p1, vertex p2, vector<vertex> g, int n, vector<vertex> &gc){
  for(int i=0; i<n; ++i){
    if( intersects(p1, p2, g[i], g[(i+1)%n]) ){
      vertex inter = intersection(p1, p2, g[i], g[(i+1)%n]);
      gc.push_back(inter);
    }
  }
}

//check if ga vertices are inside gb
void check_vertices(vector<vertex> ga, int na, vector<vertex> gb, int nb, vector<vertex> &gc, bool calc_inter){
  for(int i=0; i<na; ++i){
    if( inside(ga[i], gb, nb) )
      gc.push_back(ga[i]);
    if( calc_inter )
      check_intersect(ga[i], ga[(i+1)%na], gb, nb, gc);
  }//i
}

void sort_vertices(vector<vertex> &g, int n){
  if( n > 0 ){
    int i, j, m = 0;
    for(i=1; i<n; ++i)
      if(g[i].y < g[m].y || (g[i].y == g[m].y && g[i].x < g[m].y))
        m = i;
    swap(g[m], g[0]);
    for(i=1; i<n; ++i){
      m = i;
      for(j=i+1; j<n; ++j){
        float t = cross(g[0], g[m], g[j]);
        if( t < 0 || (t == 0 && dist(g[0], g[j]) < dist(g[0], g[m])) )
          m = j;
      }//j
      swap(g[m], g[i]);
    }//i
  }
}

float area(vector<vertex> g, int n){
  float a = 0;
  for(int i=0; i<n; ++i)
    a += (g[(i+n-1)%n].x * g[i].y) - (g[i].x * g[(i+n-1)%n].y);
  return (a < 0) ? ( a / -2 ) : (a / 2);
}

int main(){
  int na, nb, nc; //number of vertices of graph a, b, c
  while( scanf("%d", &na) && (na != 0) ){
    vector<vertex> ga, gb, gc;
    vertex temp;
    
    for(int i=0; i<na; ++i){
      scanf("%f%f", &temp.x, &temp.y);
      ga.push_back(temp);
    }
    scanf("%d", &nb);
    for(int i=0; i<nb; ++i){
      scanf("%f%f", &temp.x, &temp.y);
      gb.push_back(temp);
    }
    
    check_vertices(ga, na, gb, nb, gc, 1);
    check_vertices(gb, nb, ga, na, gc, 0);
    nc = gc.size();
    sort_vertices(gc, nc);

    printf("%8.2f", area(ga, na) + area(gb, nb) - (2 * area(gc, nc)));
  }//while
  printf("\n");
  return 0;
}