欢迎您访问程序员文章站本站旨在为大家提供分享程序员计算机编程知识!
您现在的位置是: 首页

计算几何学 | 线段相交问题 | 平面扫描 | Segment Intersections: Manhattan Geometry | C/C++实现 | 曼哈顿几何

程序员文章站 2022-04-01 13:35:11
...

问题描述

现给出n条平行于x轴或y轴的线段,请输出其交点数。

输入:
第1行输入线段数n。接下来n行输入n条线段。每条线段按照下述格式给出:
x1x_1 y1y_1 x2x_2 y2y_2
上面分别为线段两端点的坐标。各输入数据均为整数。
输出:
输出交点的总数,占1行。
限制:
1 ≤ n ≤ 100000
互相平行的2条或更多线段之间不存在重叠的点或线段
交点数不超过1000000。

输入示例

6
2 2 2 5
1 3 5 3
4 1 4 4
5 2 7 2
6 1 6 3
6 5 6 7

输出示例

3

讲解

求n条线段的交点,可以用抽选配对的方法遍历所有可能的线段对,将交点一一列举(计数)。但是复杂度为O(n2)O(n^2),不理想。

与轴平行的线段相交问题(曼哈顿几何)可以通过平面扫描(sweep)高效地求解。平面扫描算法的思路是将一条与x轴(或y轴)平行的直线向上(向右)平行移动,在移动过程中寻找交点。这条直线称为扫描线。

扫描线并不是按照固定的间隔逐行扫描,而是在每次遇到平面上线段的端点时停止移动,然后检查该位置上的线段交点。为了进行上述处理,我们需要先将输入的线段端点按照y值排序,让扫描线向y轴正方向移动。

在扫描线移动的过程中,算法会将扫描线穿过的垂直线段(与y轴平行)临时记录下来,等到扫描线与水平线段(与x轴平行)重叠时,检查水平线段的范围内是否存在垂直线段上的点,然后将这些点作为交点输出。为提高处理效率,我们可以应用二叉搜索树来保存扫描穿过的垂直线段。线段相交问题的平面扫描算法具体如下

平面扫描:
1.将已输入线段的端点按y坐标升序排列,添加至表EP
2.将二叉搜索树T置为空
3.按顺序取出EP的端点(相当于让扫描线自下而上移动),进行以下处理:
如果取出的端点为垂直线段的上端点,则从T中删除该线段的x坐标
如果取出的端点为垂直线段的下端点,则将该线段的x坐标插入T
如果取出的端点为水平线段的左端点(扫描线与水平线段重合时),将该水平线段的两端点作为搜索范围,输出T中包含的值(即垂直线段的x坐标)。

使用平衡的二叉搜索树后,1此搜索的复杂度为O(logn)O(logn),由于这个值小于2n,所以二叉树带来的复杂度为O(nlogn)O(nlogn)。算法整体的复杂度还与交点数k有关,因此本题中介绍的平面扫描算法的复杂度为O(nlogn+k)O(nlogn+k)

另外这个问题还可用区间树(segment tree)高效求解。

平面搜索:

#define BOTTOM 0
#define LEFT 1
#define RIGHT 2
#define TOP 3

class EndPoint {
	public:
		Point p;
		int seg, st;//输入线段的ID,端点的种类
		EndPoint() {}
		EndPoint(Point p, int seg, int st): p(p), seg(seg), st(st) {}
		
		bool operator < (const EndPoint &ep) const {
			//按y坐标升序排序
			if(p.y == ep.p.y) {
				return st < ep.st;//y相同时,按照下端点、左端点、右端点、上端点的顺序排列 
			} else return p.y < ep.p.y; 
		}
}; 

EndPoint EP[2 * 100000];//端点列表

//线段相交问题:曼哈顿几何
int manhattanIntersection(vector<Segment> S) {
	int n = S.size();
	
	for(int i = 0, k = 0; i < n; i++) {
		//调整端点p1、p2,保证左小右大
		if(S[i].p1.y == S[i].p2.y) {
			if(S[i].p1.x > S[i].p2.x) swap(S[i].p1, S[i].p2);
		} else if (S[i].p1.y > S[i].p2.y) swap(S[i].p1, S[i].p2);
		
		if(S[i].p1.y == S[i].p2.y) {
			EP[k++] = EndPoint(S[i].p1, i, LEFT);
			EP[k++] = EndPoint(S[i].p2, i, RIGHT);
		} else {
			EP[k++] = EndPoint(S[i].p1, i, BOTTOM);
			EP[k++] = EndPoint(S[i].p2, i, TOP);
		}
	}
	
	sort(EP, EP + (2 * n));//按端点的y坐标升序排列
	
	set<int> BT;//二叉搜索树
	BT.insert(1000000001);//设置标记
	int cnt = 0;
	
	for(int i = 0; i < 2 * n; i++) {
		if(EP[i].st == TOP) {
			BT.erase(EP[i].p.x);//删除上端点 
		} else if(EP[i].st == BOTTOM) {
			BT.insert(EP[i].p.x);//添加下端点 
		} else if(EP[i].st == LEFT) {
			set<int>::iterator b = BT.lower_bound(S[EP[i].seg].p1.x);//O(log n)
			set<int>::iterator e = BT.upper_bound(S[EP[i].seg].p2.x);//O(log n)
			cnt += distance(b, e);//加上b和e的距离(点数) O(k) 
		}
	} 
	
	return cnt;
} 

实现上述平面扫描算法时要注意各种处理的顺序,以免在一条扫描线上同时进行多个处理时遗漏交点。上述算法在一条扫描线上同时进行线段(x坐标的值)的删除、插入、搜索时,会按照下端点、左端点、右端点、上端点的顺序排列端点,从而避免这一问题。

AC代码如下

#include<iostream>
#include<cmath>
#include<vector> 
#include<algorithm>
#include<set>
using namespace std;

#define EPS (1e-10)
#define equals(a, b) (fabs((a) - (b)) < EPS)

class Point {//Point类,点 
	public:
		double x, y;
		
		Point(double x = 0, double y = 0): x(x), y(y) {}

		Point operator + (Point p) { return Point(x + p.x, y + p.y); }
		Point operator - (Point p) { return Point(x - p.x, y - p.y); }
		Point operator * (double a) { return Point(a * x, a * y); }
		Point operator / (double a) { return Point(x / a, y / a); }

		double abs() { return sqrt(norm()); }
		
		double norm() { return x * x + y * y; }
		
		bool operator < (const Point &p) const {
			return x != p.x ? x < p.x : y < p.y;
		}

		bool operator == (const Point &p) const {
			return fabs(x - p.x) < EPS && fabs(y - p.y) < EPS;
		}
};

typedef Point Vector;//Vector类,向量 

struct Segment{//Segment 线段 
	Point p1, p2;
};

double dot(Vector a, Vector b) {//内积 
	return a.x * b.x + a.y * b.y;
}

#define BOTTOM 0
#define LEFT 1
#define RIGHT 2
#define TOP 3

class EndPoint {
	public:
		Point p;
		int seg, st;//输入线段的ID,端点的种类
		EndPoint() {}
		EndPoint(Point p, int seg, int st): p(p), seg(seg), st(st) {}
		
		bool operator < (const EndPoint &ep) const {
			//按y坐标升序排序
			if(p.y == ep.p.y) {
				return st < ep.st;//y相同时,按照下端点、左端点、右端点、上端点的顺序排列 
			} else return p.y < ep.p.y; 
		}
}; 

EndPoint EP[2 * 100000];//端点列表

//线段相交问题:曼哈顿几何
int manhattanIntersection(vector<Segment> S) {
	int n = S.size();
	
	for(int i = 0, k = 0; i < n; i++) {
		//调整端点p1、p2,保证左小右大
		if(S[i].p1.y == S[i].p2.y) {
			if(S[i].p1.x > S[i].p2.x) swap(S[i].p1, S[i].p2);
		} else if (S[i].p1.y > S[i].p2.y) swap(S[i].p1, S[i].p2);
		
		if(S[i].p1.y == S[i].p2.y) {
			EP[k++] = EndPoint(S[i].p1, i, LEFT);
			EP[k++] = EndPoint(S[i].p2, i, RIGHT);
		} else {
			EP[k++] = EndPoint(S[i].p1, i, BOTTOM);
			EP[k++] = EndPoint(S[i].p2, i, TOP);
		}
	}
	
	sort(EP, EP + (2 * n));//按端点的y坐标升序排列
	
	set<int> BT;//二叉搜索树
	BT.insert(1000000001);//设置标记
	int cnt = 0;
	
	for(int i = 0; i < 2 * n; i++) {
		if(EP[i].st == TOP) {
			BT.erase(EP[i].p.x);//删除上端点 
		} else if(EP[i].st == BOTTOM) {
			BT.insert(EP[i].p.x);//添加下端点 
		} else if(EP[i].st == LEFT) {
			set<int>::iterator b = BT.lower_bound(S[EP[i].seg].p1.x);//O(log n)
			set<int>::iterator e = BT.upper_bound(S[EP[i].seg].p2.x);//O(log n)
			cnt += distance(b, e);//加上b和e的距离(点数) O(k) 
		}
	} 
	
	return cnt;
} 

int main(){
	int n;
	cin>>n;
	
	vector<Segment> S;
	Segment s;
	
	while(n--){
		cin>>s.p1.x>>s.p1.y>>s.p2.x>>s.p2.y;
		S.push_back(s);
	}
	
	cout<<manhattanIntersection(S)<<endl;
}

注:以上本文未涉及代码的详细解释参见:计算几何学