GDAL出现 TopologyException: Input geom 1 is invalid: Self-intersection at or near point…如何解决
最近一段时间在写关于GDAL的代码,我在对两个Polygon求交时出现这个问题:
百度后得知是多边形内部产生了包含其中的一个小多边形导致,给出的方法是求距离为0的缓冲区,但使用后问题没有得到解决。
于是,我输出了多边形全部的顶点(输出多边形顶点网上有代码,这里就不粘了),根据代码报错给的提示找出现问题的顶点(719273.8603311308, 2543999.7224521482)。
发现原来这个顶点重复了三次!那么原则上只要删除两个,留一个即可解决问题,于是我写了一个算法:
//判断OGRPoint指针指向的点是否相同
bool equalOGRPoint(OGRPoint *oLeftPoint, OGRPoint *oRightPoint)
{
if (fabs(oLeftPoint->getX() - oRightPoint->getX()) < 10e-5)
{
if (fabs(oLeftPoint->getY() - oRightPoint->getY()) < 10e-5)
{
return true;
}
}
return false;
}
//删除错误顶点,生成新多边形返回
OGRGeometry * dealWithSelfIntersection(OGRGeometry * geom)
{
//重复点的位置
int flag = 0;
OGRPoint *oPoint = new OGRPoint();
OGRPoint *oPointNext = new OGRPoint();
OGRLinearRing *oLineString = ((OGRPolygon *)geom)->getExteriorRing();
//注意闭合环首位顶点坐标是相同的,五个顶点的多边形实则有六个顶点坐标
for (int i = 0; i < oLineString->getNumPoints()-3; i++)
{
oLineString->getPoint(i, oPoint);
oLineString->getPoint(i + 2, oPointNext);
if (equalOGRPoint(oPoint, oPointNext))
{
flag = i+1;
break;
}
}
//GDAL没有删点的函数,创建一个新的多边形返回
OGRPolygon* pPolygon = (OGRPolygon*)OGRGeometryFactory::createGeometry(wkbPolygon);
//新多边形边界
OGRLinearRing *oLineStringNew = new OGRLinearRing();
for (int i = 0; i < oLineString->getNumPoints() - 1; i++)
{
if (i != (flag + 1) && i != (flag + 2))
{
oLineString->getPoint(i, oPoint);
oLineStringNew->addPoint(oPoint);
}
}
oLineStringNew->closeRings();
pPolygon->addRing(oLineStringNew);
return pPolygon;
}
调用很简单,再求交前,对求交的OGRGeometry* 调用该函数
if (pPolygon == NULL)
{
pPolygon= dealWithSelfIntersection(pPolygon);
}
即可
当然,在我遇到的这个问题里,多边形并不是生成了一个环,而是顶点重复三次,出现了类似内环的效果。不过上述代码只要是三个顶点的内环都是可以使用的。如果内环顶点出现多个的情况,代码略微更改即可!
如果大家遇到与GDAL(C++)有关的问题,欢迎给我发邮件 Chengtcug@foxmail.com !