Weiler Atherton 算法

Weiler Atherton 算法是通用的多边形裁减算法，可以用来求两个多边形的交集。在使用kitti_eval时评价AOS的部分，是借助boost.geometry求两个bounding box的交， 它的实现就是基于Weiler Atherton 算法。

算法步骤：
1、首先需要准备两个多边形顶点集，这些点可以是顺时针或逆时针排列
2、被裁减的主多边形记为P，裁减窗口多边形记为Q，求出P，Q所有的交点，并按序插入到原有的顶点集中


3、将交点分为两类：
　　一类称入点，即被裁剪多边形由此点进入裁剪窗口
　　一类称出点，即被裁剪多边形由此点离开裁剪窗口
　可以看出，相交区域是由PQ出入点交替组成(中间可以由原顶点)。

4、交替遍历即可得到相交区域

代码实现：

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

#define Size 600
using namespace std;

typedef float Color[3];
struct Point 
{
	int x, y;
	Point() = default;
	Point(int x1, int y1) : x(x1), y(y1){}
};
typedef struct IntersectionPoint 
{
	int pointFlag;
	int index0, index1;
	Point p;
	bool inFlag;
	int dis;
}IP;

class Pg 
{
public:
	vector pts;
	Pg(void);
	~Pg(void);
	void drawPgLine(Color c);
};
Pg::Pg(void) 
{
}
Pg::~Pg(void) 
{
}
void Pg::drawPgLine(Color c) 
{
	glColor3fv(c);
	glLineWidth(2.0);
	glBegin(GL_LINE_LOOP);
	int size = pts.size();
	for (int i = 0; i < size; i++)
		glVertex2i(pts[i].x, pts[i].y);
	glEnd();
}
// 判断是否在多边形内
bool isPointInsidePg(Point p, Pg& py) 
{
	int cnt = 0, size = py.pts.size();
	for (int i = 0; i < size; i++) {
		Point p1 = py.pts[i];
		Point p2 = py.pts[(i + 1) % size];
		if (p1.y == p2.y) continue;
		if (p.y < min(p1.y, p2.y)) continue;
		if (p.y >= max(p1.y, p2.y)) continue;
		double x = (double)(p.y - p1.y) * (double)(p2.x - p1.x) / (double)(p2.y - p1.y) + p1.x;
		if (x > p.x) cnt++;
	}
	return (cnt % 2 == 1); 
}

int cross(Point& p0, Point& p1, Point& p2) 
{
	return ((p2.x - p0.x) * (p1.y - p0.y) - (p1.x - p0.x) * (p2.y - p0.y));
}
bool onSegment(Point& p0, Point& p1, Point& p2) 
{
	int minx = min(p0.x, p1.x), maxx = max(p0.x, p1.x);
	int miny = min(p0.y, p1.y), maxy = max(p0.y, p1.y);
	if (p2.x >= minx && p2.x <= maxx && p2.y >= miny && p2.y <= maxy) return true;
	return false;
}

bool segmentsIntersect(Point& p1, Point& p2, Point& p3, Point& p4) 
{
	int d1 = cross(p3, p4, p1);
	int d2 = cross(p3, p4, p2);
	int d3 = cross(p1, p2, p3);
	int d4 = cross(p1, p2, p4);
	if (((d1 > 0 && d2 < 0) || (d1 < 0 && d2 > 0)) &&
		((d3 > 0 && d4 < 0) || (d3 < 0 && d4 > 0)))
		return true;
	if (d1 == 0 && onSegment(p3, p4, p1)) return true;
	if (d2 == 0 && onSegment(p3, p4, p2)) return true;
	if (d3 == 0 && onSegment(p1, p2, p3)) return true;
	if (d4 == 0 && onSegment(p1, p2, p4)) return true;
	return false;
}
Point getintersectPoint(Point p1, Point p2, Point p3, Point p4) 
{
	Point p;
	int b1 = (p2.y - p1.y) * p1.x + (p1.x - p2.x) * p1.y;
	int b2 = (p4.y - p3.y) * p3.x + (p3.x - p4.x) * p3.y;
	int D = (p2.x - p1.x) * (p4.y - p3.y) - (p4.x - p3.x) * (p2.y - p1.y);
	int D1 = b2 * (p2.x - p1.x) - b1 * (p4.x - p3.x);
	int D2 = b2 * (p2.y - p1.y) - b1 * (p4.y - p3.y);
	p.x = D1 / D;
	p.y = D2 / D;
	return p;
}
// 获得两组多边形的交点
void generateIntersectPoints(Pg& pyclip, Pg& py, list& iplist) 
{
	int clipSize = pyclip.pts.size(), pySize = py.pts.size();

	for (int i = 0; i < clipSize; i++) {
		Point p1 = pyclip.pts[i];
		Point p2 = pyclip.pts[(i + 1) % clipSize]; 
		for (int j = 0; j < pySize; j++) {
			Point p3 = py.pts[j];
			Point p4 = py.pts[(j + 1) % pySize];
			if (segmentsIntersect(p1, p2, p3, p4)) {
				IP ip;
				ip.index0 = j;
				ip.index1 = i;
				ip.p = getintersectPoint(p1, p2, p3, p4);
				iplist.push_back(ip); 
			}
		}
	}
}
int getDistance(Point& p1, Point& p2) {
	return (p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y);
}
bool distanceComparator(IP& ip1, IP& ip2) {
	return ip1.dis < ip2.dis;
}
void generateList(Pg& py, list& iplist, list& comlist, int index) {
	int size = py.pts.size(); // 原来的点数
	list::iterator it; // 

	for (int i = 0; i < size; i++) {
		Point p1 = py.pts[i];
		IP ip;
		ip.pointFlag = 0;
		ip.p = p1;
		comlist.push_back(ip);// 
		list oneSeg;
		// 遍历iplist 交点
		for (it = iplist.begin(); it != iplist.end(); it++) {
			if ((index == 0 && i == it->index0) ||
				(index == 1 && i == it->index1)) {
				it->dis = getDistance(it->p, p1);// 
				it->pointFlag = 1;
				oneSeg.push_back(*it); // 
			}
		}
		oneSeg.sort(distanceComparator);// 按序排列，可能存在多个交点
		for (it = oneSeg.begin(); it != oneSeg.end(); it++) // 
			comlist.push_back(*it);
	}
}
void getPgPointInOut(list& Pglist, Pg& pyclip) {
	bool inFlag;
	list::iterator it;
	for (it = Pglist.begin(); it != Pglist.end(); it++) {
		if (it->pointFlag == 0) { // 
			if (isPointInsidePg(it->p, pyclip))
				inFlag = true;
			else inFlag = false;
		}
		else {
			inFlag = !inFlag;
			it->inFlag = inFlag;//
		}
	}
}
bool operator==(Point& p1, Point& p2) {
	return p1.x == p2.x && p1.y == p2.y;
}
// 为交点赋上出入点属性
void getClipPointInOut(list& cliplist, list& Pglist) {
	list::iterator it, it1;
	for (it = cliplist.begin(); it != cliplist.end(); it++) {
		if (it->pointFlag == 0) continue;
		for (it1 = Pglist.begin(); it1 != Pglist.end(); it1++) {
			if (it1->pointFlag == 0) continue; // 原顶点
			if (it->p == it1->p) it->inFlag = it1->inFlag;
		}
	}
}
// 计算多边形面积
double calculate_ploygon_area(const vector& points){
	//points 应当是顺时针或逆时针排列
	int n = points.size();
	double sum1 = 0.0, sum2 = 0.0;
	for(uint i = 0;i& Pglist, list& cliplist) {
	double intersection_area = 0.0;
	list::iterator it, it1;
	Pg py;
	Color c = { 0.0, 0.0, 1.0 };

	for (it = Pglist.begin(); it != Pglist.end(); it++)
		if (it->pointFlag == 1 && it->inFlag) break; // 找到第一个inFlag=true的交点
	py.pts.clear();
	auto it_begin = it;
	while (true) {
		//py.pts.push_back(it->p);
		for (; ; it++) {
			if(it==Pglist.end()) it = Pglist.begin();
			if (it->pointFlag == 1 && !it->inFlag) break;
			py.pts.push_back(it->p);
		}
		for (it1 = cliplist.begin(); it1 != cliplist.end(); it1++){
			if(it1==cliplist.end()) it1 = cliplist.begin();
			if (it1->p == it->p) break; // it1追到it
		}
		for (; ; it1++) {
			if(it1==cliplist.end()) it1 = cliplist.begin();
			if (it1->pointFlag == 1 && it1->inFlag) break;
			py.pts.push_back(it1->p);
		}

		if (py.pts[0]==it1->p) {
			intersection_area += calculate_ploygon_area(py.pts);
			py.drawPgLine(c);
			py.pts.clear();
			for (; it != Pglist.end(); it++)
				if (it->pointFlag == 1 && it->inFlag) break; //寻找下一个环

			if (it == it_begin||it==Pglist.end()) break;// 退出条件
			else continue;
		}else{
			for (; ; it++)
			{
				if(it==Pglist.end()) it = Pglist.begin();
				if (it->p == it1->p) break; // it 追到it1
			}
				
		}

	}
	std::cout<<"the intersection area is: " < iplist, Pglist, cliplist;
	generateIntersectPoints(pyclip, py, iplist); // 获取两个多边形的相交点
	generateList(py, iplist, Pglist, 0);//

	generateList(pyclip, iplist, cliplist, 1); // 

	getPgPointInOut(Pglist, pyclip); // 
	getClipPointInOut(cliplist, Pglist); // 将内外点属性赋给cliplist

	generateClipArea(Pglist, cliplist);
}

void init()
{
	glClearColor(0.0, 0.0, 0.0, 0.0);
	glColor3f(1.0, 0.0, 0.0);
	glPointSize(1.0);
	glMatrixMode(GL_PROJECTION);
	glLoadIdentity();
	gluOrtho2D(0.0, Size - 1, 0.0, Size - 1);
}

void GenerateRandomSimplePg(Pg &G, int M)
{
	Point P;
	G.pts.clear();
	for (int i = 0; i < M; ++i)
	{
		bool flag;
		do
		{
			P.x = rand() % Size;
			P.y = rand() % Size;
			flag = true;
			for (int j = 1; j < i - 1; ++j)
				if (segmentsIntersect(G.pts[j - 1], G.pts[j], G.pts[i - 1], P))
				{
					flag = false;
					break;
				}
			if (flag && i == M - 1)
			{
				for (int j = 2; j < i; ++j)
					if (segmentsIntersect(G.pts[j - 1], G.pts[j], P, G.pts[0]))
					{
						flag = false;
						break;
					}
			}
		} while (!flag);
		G.pts.push_back(P);
	}
}

void KeyboardAction(unsigned char key, int x, int y)
{
	if(key=='q')
		exit(0);
}

void display()
{
	glClear(GL_COLOR_BUFFER_BIT);
	glEnable(GL_POINT_SMOOTH);

	Pg pyclip, py;

	GenerateRandomSimplePg(pyclip, 4);
	GenerateRandomSimplePg(py, 4);

	// Point p1, p2, p3, p4; 
	// p1.x = 100, p1.y = 100;
	// p2.x = 500, p2.y = 100;
	// p3.x = 500, p3.y = 500;
	// p4.x = 100, p4.y = 500;
	// pyclip.pts.push_back(p1);
	// pyclip.pts.push_back(p2);
	// pyclip.pts.push_back(p3);
	// pyclip.pts.push_back(p4);

	// Point p5, p6, p7, p8;
	// p5.x = 150, p5.y = 252;
	// p6.x = 400, p6.y = 53;
	// p7.x = 550, p7.y = 510;
	// //p8.x = 68, p8.y = 245;
	// py.pts.push_back(p5);
	// py.pts.push_back(p6);
	// py.pts.push_back(p7);
	// //py.pts.push_back(p8);

	int size = pyclip.pts.size();
	for (int i = 0; i < size; ++i)
		cout << pyclip.pts[i].x << " " << pyclip.pts[i].y << endl;
	cout << endl;
	size = py.pts.size();
	for (int i = 0; i < size; ++i)
		cout << py.pts[i].x << " " << py.pts[i].y << endl;

	Color a = { 1.0, 0.0, 0.0 };
	Color b = { 0.0, 1.0, 0.0 };
	py.drawPgLine(a);
	pyclip.drawPgLine(b);
	weilerAtherton(pyclip, py); // 算法入口

	glFlush();
}

int main(int argc, char **argv)
{
	srand(time(NULL));

	glutInit(&argc, argv);
	glutInitDisplayMode(GLUT_SINGLE | GLUT_RGB);
	glutInitWindowSize(Size, Size);
	glutInitWindowPosition(100, 100);
	glutCreateWindow("Weiler-Atherton Clipping Algorithm");

	glutKeyboardFunc(KeyboardAction);
	glutDisplayFunc(display); // 函数入口

	init();

	glutMainLoop();

	return 0;
}

编译：g++ WA-clipping.cpp -o WA -lGL -lGLU -lglut
// 这个代码现在可能还有点bug的，使用随机数有时不能显示

参考：
https://github.com/Espade/Weiler-Atherton-Clipping