两个多边形求交的实现需要几个模块 (cgal中有insect函数,但是必须要求使用CGAL::Exact_predicates_exact_constructions_kernel) 本人出于其他想法,没有把基于Exact_predicates_inexact_constructions_kernel核的Polygon_2转换为 Exact_predicates_exact_constructions_kernel的核,而是自己写了一个求交函数自用.
1.线线求交 1.1 cgal版
代码语言:javascript复制//cgal实现版本
bool GetInsertionOfLineAndLine(double bsx,double bsy,double bex,double bey, double sx, double sy, double ex, double ey, double&rx, double &ry) {
CGAL::Exact_predicates_exact_constructions_kernel::Segment_2 lin(CGAL::Exact_predicates_exact_constructions_kernel::Point_2(bsx, bsy), CGAL::Exact_predicates_exact_constructions_kernel::Point_2(bex, bey));
CGAL::Exact_predicates_exact_constructions_kernel::Segment_2 seg(CGAL::Exact_predicates_exact_constructions_kernel::Point_2(sx, sy), CGAL::Exact_predicates_exact_constructions_kernel::Point_2(ex, ey));
CGAL::cpp11::result_of<CGAL::Exact_predicates_exact_constructions_kernel::Intersect_2(CGAL::Exact_predicates_exact_constructions_kernel::Segment_2, CGAL::Exact_predicates_exact_constructions_kernel::Segment_2)>::type
result = intersection(lin, seg);
if (result) {
if (const CGAL::Exact_predicates_exact_constructions_kernel::Segment_2* s = boost::get<CGAL::Exact_predicates_exact_constructions_kernel::Segment_2>(&*result)) {
}
else {
const CGAL::Exact_predicates_exact_constructions_kernel::Point_2* p = boost::get<CGAL::Exact_predicates_exact_constructions_kernel::Point_2 >(&*result);
rx = to_double((*p).x());
ry = to_double((*p).y());
//std::cout << jx << ";" << jy << endl;
return true;
}
}
return false;
}
2.判断点在多边形内 2.1通用版(判断内角和是否等于360°)
代码语言:javascript复制double GetAngle(Point_2 p0, Point_2 p1, Point_2 p2)
{
double Pi = 3.1415926;
double EPS = 0.000001;
Point_2 a(p0.x() - p1.x(), p0.y() - p1.y());
Point_2 b(p2.x() - p1.x(), p2.y() - p1.y());
double angle1 = atan2(a.y(), a.x()) * 180 / Pi;
double angle2 = atan2(b.y(), b.x()) * 180 / Pi;
if (angle1*angle2 > -EPS)
return abs(angle1 - angle2)*Pi / 180;
else {
double tmp = abs(angle1) abs(angle2);
if (tmp - 180 > EPS)
tmp = 360 - tmp;
return tmp * Pi / 180;
}
}
bool IsPointInPolygon(const Point_2 p, const Polygon_2 poly)
{
double EPS = 0.000001;
double Pi = 3.1415926;
int length = poly.size();
double angle = 0;
for (int i = 0; i < length; i)
{
angle = GetAngle(poly.vertex(i), p, poly.vertex((i 1) % length));
}
return abs(angle - 2 * Pi) < EPS;
}
bool IsPolygonInPolygon(const Polygon_2 p, const Polygon_2 poly)
{
int len = p.size();
for (int i = 0; i < len; i)
{
if (!IsPointInPolygon(p.vertex(i), poly))
return false;
}
return true;
}
3.多边形求交 3.1 通用版(Polygon逆时针存储vertex) a.通过线线求交找到两个交点 p1,p2 以及对应的相交边的起始顶点 i1,j1 和 i2,j2 b.通过判断 i1 1是否在polyb中,获取相交部分的逆时针索引vector
代码语言:javascript复制bool Intersect(const Polygon_2 polya, const Polygon_2 polyb, Polygon_2 *intersect)
{
int lena = polya.size();
int lenb = polyb.size();
double *_polyax = new double[lena];
double *_polyay = new double[lena];
double *_polybx = new double[lenb];
double *_polyby = new double[lenb];
for (int i = 0; i < lena; i)
{
_polyax[i] = polya.vertex(i).x();
_polyay[i] = polya.vertex(i).y();
}
for (int i = 0; i < lenb; i)
{
_polybx[i] = polyb.vertex(i).x();
_polyby[i] = polyb.vertex(i).y();
}
bool first = false;
bool second = false;
int i1, i2, j1, j2;
double p1x, p1y, p2x, p2y;
bool IsInsertion = false;
for (int i = 0; i < lena; i)
{
for (int j = 0; j < lenb; j)
{
IsInsertion = false;
double jx, jy;
IsInsertion = GetInsertionOfLineAndLine(_polyax[i], _polyay[i], _polyax[(i 1)%lena], _polyay[(i 1) % lena], _polybx[j], _polyby[j], _polybx[(j 1) % lenb], _polyby[(j 1) % lenb], jx, jy);
if (IsInsertion) {
if (!first) {
first = true;
i1 = i; j1 = j;
p1x = jx; p1y = jy;
}
else {
second = true;
i2 = i; j2 = j;
p2x = jx; p2y = jy;
}
break;
}
}
}
if (!first || !second)
return false;
Polygon_2 result; int s1, s2;
//if vertex(i1 1) is located in poly2,get p1->i1 1->...->i2->p2->j2 1->...->j1
if (this->IsPointInPolygonByBao(polya.vertex((i1 1)%lena), polyb))
{
s1 = (i1 1) % lena;
s2 = (j2 1) % lenb;
result.push_back(Point_2(p1x, p1y));
if (s1 <= i2)
{
for (int i = s1; i <= i2; i)
result.push_back(polya.vertex(i));
}
else {
for (int i = s1; i < lena; i)
result.push_back(polya.vertex(i));
for(int i=0;i<=i2; i)
result.push_back(polya.vertex(i));
}
result.push_back(Point_2(p2x, p2y));
if(j1 != j2)
if (s2 <= j1)
{
for (int i = s2; i <= j1; i)
result.push_back(polyb.vertex(i));
}
else {
for (int i = s2; i < lenb; i)
result.push_back(polyb.vertex(i));
for (int i = 0; i <= j1; i)
result.push_back(polyb.vertex(i));
}
}
//if vertex(i1 1) isn't located in poly2,get p1->j1 1->...->j2->p2->i2 1->...->i1
else {
s1 = (j1 1) % lenb;
s2 = (i2 1) % lena;
result.push_back(Point_2(p1x, p1y));
if(j1!=j2)
if (s1 <= j2)
{
for (int i = s1; i <= j2; i)
result.push_back(polyb.vertex(i));
}
else {
for (int i = s1; i < lenb; i)
result.push_back(polyb.vertex(i));
for (int i = 0; i <= j2; i)
result.push_back(polyb.vertex(i));
}
result.push_back(Point_2(p2x, p2y));
if (s2 <= i1)
{
for (int i = s2; i <= i1; i)
result.push_back(polya.vertex(i));
}
else {
for (int i = s2; i < lena; i)
result.push_back(polya.vertex(i));
for (int i = 0; i <= i1; i)
result.push_back(polya.vertex(i));
}
}
*intersect = result;
return true;
}