Fast algorithm for union of lwpolyline-a and lwpolyline-b

[1127204185]

Fast algorithm for union of lwpolyline-a and lwpolyline-b

[tbrammer · Accepted solution]

You could create two AcDbRegions region1 and region2 from the polylines and unite them with

   region1->booleanOper(AcDb::kBoolUnite, region2);

Then analyze the resulting region1 using the BREP library (<ARX>\utils\brep) which allows you to trace inner and outer boundaries.

Try to build the <ARX>\utils\brep\samples\brepsamp\brsample.vcxproj, load the brsample.arx and enter the command BRDUMP. Hit [Return] until you are prompted to "Pick a solid" and select the region.

The command will dump the topology of the region organized like this:

Complex ( 1 Complex for each separated area)

Shell   (in a region each Complex has 1 Shell)

Face    (area with 0..N holes that has no common edge with another area. Faces may touch at single points)

Loop   (outer or inner border of a face. An "interior loop" is a hole. An "exterior loop" is the outer border)

Vertex (a point)

If you have "touching faces" like in your picture, you might probably get two separate faces. In this case you have to look for vertices that are joint by two faces. You can trace the outline by following loops but switching to a different face when you find a joint vertex.

I have attached a DWG with two regions and the output of BRDUMP for the two regions.

[1127204185]

ads_name ss;
struct resbuf *rb=acutBuildList(RTDXF0,_T("REGION"),RTNONE);
acutPrintf(_T("\n请选择两个需要布尔运算(交集)的面域: "));
if(RTNORM != acedSSGet(NULL,NULL,NULL,rb,ss))
{
acutRelRb(rb);
return;
}
acutRelRb(rb);
Adesk::Int32 nssLen=0;
acedSSLength(ss,&nssLen);
if(nssLen<2)
{
acedSSFree(ss);
return;
}
ads_name ent;
acedSSName(ss,0,ent);
AcDbObjectId objId1;
acdbGetObjectId(objId1,ent);
acedSSName(ss,1,ent);
AcDbObjectId objId2;
acdbGetObjectId(objId2,ent);
acedSSFree(ss);
//使用ARX智能指针写打开面域对象
AcDbObjectPointer<AcDbRegion>pRegion1(objId1,AcDb::kForWrite);
if (Acad::eOk != pRegion1.openStatus())
{
return;
}
AcDbObjectPointer<AcDbRegion>pRegion2(objId2,AcDb::kForWrite);
if (Acad::eOk != pRegion2.openStatus())
{
return;
}
//执行布尔操作
Acad::ErrorStatus es= pRegion1->booleanOper(AcDb::kBoolIntersect,pRegion2);
if (Acad::eOk != eOk)
{
acutPrintf(_T("\n布尔操作失败!,错误码:%s"),acadErrorStatusText(es));
return;
}
//判断是否有交集，如果为空,则删除对象。
if (pRegion1->isNull())
{
acutPrintf(_T("\n没有交集,删除对象!"));
pRegion1->erase();
pRegion2->erase();
return;
}
//删除多余的面域实体，执行交集布尔运算后，参数计算的面域面积为0。
if (pRegion2->isNull())
{
pRegion2->erase();
}

//打开当前数据库的当前块表记录
AcDbBlockTableRecordPointer pBlkRcd(curDoc()->database()->currentSpaceId(),AcDb::kForWrite);
if (Acad::eOk != pBlkRcd.openStatus())
{
acutPrintf(_T("\n打开块表记录失败!,错误码:%s"),acadErrorStatusText(pBlkRcd.openStatus()));
return;
}
//分解面域
AcDbVoidPtrArray ptrArray;
es = pRegion1->explode(ptrArray);
if(Acad::eOk != es)
{
acutPrintf(_T("\n面域分解操作失败!,错误码:%s"),acadErrorStatusText(es));
return;
}
//遍历添加分解后的实体到当前空间,如果不加入数据库就需要自己delete分解后的对象
for (int i=0;i<ptrArray.length();i++)
{
AcDbEntity*pEnt=(AcDbEntity*)ptrArray.at(i);
if (pEnt!=NULL)
{
if (pEnt->isKindOf(AcDbRegion::desc()))
{
AcDbVoidPtrArray subptrArray;
es = pEnt->explode(subptrArray);
if(Acad::eOk != es)
{
acutPrintf(_T("\n子面域分解操作失败!,错误码:%s"),acadErrorStatusText(es));
continue;
}
for (int j=0;j<subptrArray.length();j++)
{
AcDbEntity*pSubEnt=(AcDbEntity*)subptrArray.at(j);
if (pSubEnt!=NULL)
{
pBlkRcd->appendAcDbEntity(pSubEnt);
//关闭pSubEnt对象
pSubEnt->close();
}
}
}
else
{
pBlkRcd->appendAcDbEntity(pEnt);
}
//关闭pEnt对象
pEnt->close();
}
}
//删除面域
//pRegion1->erase();
/*
//获取拉伸点集合，(另类的获取面域端点的方式)
AcGePoint3dArray pts;
pRegion1->getStretchPoints(pts);
if (pts.length()>0)
{
for (int i=0;i<pts.length();i++)
{
AcGePoint3d pt=pts.at(i);
acutPrintf(_T("\n x=%lf,y=%lf,z=%lf"),pt.x,pt.y,pt.z);
}
}
*/

[1127204185]

I want to find another way, get the set of points, and then find the boundary

[tbrammer]

If you are looking for a more general algorithm you can do it like this:

Let's call the polylines A and B.

1. Make sure that A and B run counter clockwise (ccw). Reverse plines that run cw.
2. Calculate all section points P[0]..P[N] of A and B.
3. For each point in P[] calculate the curve paramters of p on A and B:
   A.getParamAtPoint(P[i], pA[i])
   B.getParamAtPoint(P[i], pB[i])
   Now we are going to trace the loops. Begin with i0=0.
4. Start at P[i0] on A with parameter pA[i0]
5. Follow A to the next section point P[i1] that "turn right" to B with the smalles pA[i1] having  pA[i1] > pA[i0]. If there is no pA[i1] > pA[i0] use the smalles pA[i1]
6. Switch over to B . Here you are at parameter pB[i1].
7. Like in 5. now follow B to the next section point P[i2] that turns right to A with pB[i2] > pB[i1]
8. Switch over to A. Here you are at parameter pA[i2].
9. Repeat steps 5.-8. until you return to a point with index ik you "visited" before.
   You found a closed loop P[ik], P[ik+1], P[ik+2],..., P[ik].
   This will be an inner loop (hole) if its path runs cw.  If its path runs ccw it is the exteriour shape.
10. If you didn't "visit" all N points in the loop 5.-8. start again with any "unvisited" value for i0.
    Repeat until all points are visited or until you found the exterior loop.

For an algorithm to check the orientation of a polgon see here.

To determine whether switching from A->B is a "turn right" at P[i] calculate the curve directions of A and B at P[i] by  A->getFirstDeriv(pA[i], va); B->getFirstDeriv(pB[i], vb); 

if (va.crossProduct(vb).z < 0) A->B would be a "turn right".

We are looking for right turns because they will extent an area with ccw boundaries while left turns would make them smaller.

[1127204185]

Do you want to find the intersection？

[1127204185]

Do you want to find the intersection？

[1127204185]

[1127204185]