- 积分
- 24557
- 明经币
- 个
- 注册时间
- 2004-3-17
- 在线时间
- 小时
- 威望
-
- 金钱
- 个
- 贡献
-
- 激情
-
|
楼主 |
发表于 2010-7-5 14:38:00
|
显示全部楼层
四、最小包围圆
相关的讨论链接:
http://bbs.mjtd.com/forum.php?mod=viewthread&tid=55997
http://www.objectarx.net/forum.php?mod=viewthread&tid=4870
相关的概念和算法描述就直接贴过来了,:)
以下部分摘录自:汪卫,王文平,汪嘉业,求一个包含点集所有点的最小圆的算法,软件学报,2000,11(9):1237-1240
求一个最小圆包含给定点集所有点的问题是人们在实践和理论上都十分感兴趣的一个问题.由于这个圆的圆心是到点集中最远点最近的一个点,因而在规划某些设施时很有实用价值.
这个圆心也可看成是点集的中心.此外,在图形学中,圆也常可取作边界盒,使用它可减少很多不必要的计算.
在空间数据库中可将该问题用于建立空间数据的索引以提高查询速度[1,2].这个问题看起来十分简单,但用直观的算法去解此问题,其复杂性可达O(n4),其中n为点集中点的数目.
有关此问题的讨论在计算几何的专著及论文中未见报道[3~5].
本文提出了一种新的算法并证明了这种算法的时间复杂性为O(|lg(d/R)|*n),其中R是所求的最小圆的半径,d为点集中不在圆周上但距圆周最近的点到圆周的距离.
1、算 法
第1步.在点集中任取3点A,B,C.
第2步.作一个包含A,B,C三点的最小圆.圆周可能通过这3点(如图1所示),也可能只通过其中两点,但包含第3点.后一种情况圆周上的两点一定是位于圆的一条直径的两端.
第3步.在点集中找出距离第2步所建圆圆心最远的点D.若D点已在圆内或圆周上,则该圆即为所求的圆,算法结束.否则,执行第4步.
第4步.在A,B,C,D中选3个点,使由它们生成的一个包含这4点的圆为最小.这3点成为新的A,B和C,返回执行第2步. 若在第4步生成的圆的圆周只通过A,B,C,D中的两点,则圆周上的两点取成新的A和B,从另两点中任取一点作为新的C.
2、算法正确性
本节要证明上述算法一定能终止,且最后一次求得的圆即是所要求的包含点集所有点的最小圆.
引理1.算法第4步所生成的圆的半径随着迭代过程递增. 证明:因为第4步每一次生成的圆是包含原来的A,B,C三点,又要加上圆外的一点,而上一次生成的圆是包含A,B,C的最小圆,因而新圆的半径一定比原来的圆半径要大.
定理1.上述算法是收敛的,且最后得到包含点集所有点的最小圆. 证明:因为在点集中任取3点或两点生成的包含这3点或两点的最小圆的个数是有限的.由引理1可知,算法进行过程中所生成的圆的半径是递增的,因而经过有限次迭代后,可求得这些最小圆中半径最大的一个.从算法第3步可知,只有当点集中所有点都在由3个点或两个点生成的最小圆内时,算法才结束,因而最后得到的圆一定是包含点集中所有点的最小圆.
代码-
- using System;
- using System.Collections.Generic;
- using System.Linq;
- using System.Text;
- using TlsCad.Collections;
- using Autodesk.AutoCAD.Geometry;
- namespace TlsCad.ExtendMethods
- {
- public static class PointEx
- {
- #region PointList
- /// <summary>
- /// 判断点是否属于圆
- /// </summary>
- /// <param name="ca2d"></param>
- /// <param name="pnt"></param>
- /// <returns></returns>
- public static bool IsPartOf(this CircularArc2d ca2d, Point2d pnt)
- {
- return ca2d.IsOn(pnt) || ca2d.IsInside(pnt);
- }
- /// <summary>
- /// 按两点返回最小包围圆
- /// </summary>
- /// <param name="pt1"></param>
- /// <param name="pt2"></param>
- /// <param name="ptlst">输出圆上的点</param>
- /// <returns></returns>
- public static CircularArc2d GetMinCircle(Point2d pt1, Point2d pt2, out LoopList<Point2d> ptlst)
- {
- ptlst = new LoopList<Point2d> { pt1, pt2 };
- return
- new CircularArc2d
- (
- (pt1 + pt2.GetAsVector()) / 2,
- pt1.GetDistanceTo(pt2) / 2
- );
- }
- /// <summary>
- /// 按三点返回最小包围圆
- /// </summary>
- /// <param name="pt1"></param>
- /// <param name="pt2"></param>
- /// <param name="pt3"></param>
- /// <param name="ptlst">输出圆上的点</param>
- /// <returns></returns>
- public static CircularArc2d GetMinCircle(Point2d pt1, Point2d pt2, Point2d pt3, out LoopList<Point2d> ptlst)
- {
- ptlst =
- new LoopList<Point2d> { pt1, pt2, pt3 };
- //遍历各点与下一点的向量长度,找到距离最大的两个点
- double maxLength;
- LoopListNode<Point2d> maxNode =
- ptlst.GetNodes().FindByMaxKey
- (
- out maxLength,
- node => node.Value.GetDistanceTo(node.Next.Value)
- );
- //以两点做最小包围圆
- LoopList<Point2d> tptlst;
- CircularArc2d ca2d =
- GetMinCircle(maxNode.Value, maxNode.Next.Value, out tptlst);
- //如果另一点属于该圆
- if (ca2d.IsPartOf(maxNode.Previous.Value))
- {
- //返回
- ptlst = tptlst;
- return ca2d;
- }
- //否则按三点做圆
- ptlst.SetFirst(maxNode);
- ca2d = new CircularArc2d(pt1, pt2, pt3);
- ca2d.SetAngles(0, Math.PI * 2);
- return ca2d;
- }
- /// <summary>
- /// 按四点返回最小包围圆
- /// </summary>
- /// <param name="pt1"></param>
- /// <param name="pt2"></param>
- /// <param name="pt3"></param>
- /// <param name="pt4"></param>
- /// <param name="ptlst">输出圆上的点</param>
- /// <returns></returns>
- public static CircularArc2d GetMinCircle(Point2d pt1, Point2d pt2, Point2d pt3, Point2d pt4, out LoopList<Point2d> ptlst)
- {
- LoopList<Point2d> iniptlst =
- new LoopList<Point2d> { pt1, pt2, pt3, pt4 };
- ptlst = null;
- CircularArc2d ca2d = null;
- //遍历C43的组合,环链表的优势在这里
- foreach (LoopListNode<Point2d> firstNode in iniptlst.GetNodes())
- {
- //获取各组合下三点的最小包围圆
- LoopListNode<Point2d> secondNode = firstNode.Next;
- LoopListNode<Point2d> thirdNode = secondNode.Next;
- LoopList<Point2d> tptlst;
- CircularArc2d tca2d = GetMinCircle(firstNode.Value, secondNode.Value, thirdNode.Value, out tptlst);
- //如果另一点属于该圆,并且半径小于当前值就把它做为候选解
- if (tca2d.IsPartOf(firstNode.Previous.Value))
- {
- if (ca2d == null || tca2d.Radius < ca2d.Radius)
- {
- ca2d = tca2d;
- ptlst = tptlst;
- }
- }
- }
- //返回直径最小的圆
- return ca2d;
- }
- /// <summary>
- /// 按点集返回最小包围圆
- /// </summary>
- /// <param name="pnts"></param>
- /// <param name="ptlst">输出圆上的点</param>
- /// <returns></returns>
- public static CircularArc2d GetMinCircle(this List<Point2d> pnts, out LoopList<Point2d> ptlst)
- {
- //点数较小时直接返回
- switch (pnts.Count)
- {
- case 0:
- ptlst = new LoopList<Point2d>();
- return null;
- case 1:
- ptlst = new LoopList<Point2d> { pnts[0] };
- return new CircularArc2d(pnts[0], 0);
- case 2:
- return GetMinCircle(pnts[0], pnts[1], out ptlst);
- case 3:
- return GetMinCircle(pnts[0], pnts[1], pnts[2], out ptlst);
- case 4:
- return GetMinCircle(pnts[0], pnts[1], pnts[2], pnts[3], out ptlst);
- }
- //按前三点计算最小包围圆
- Point2d[] tpnts = new Point2d[4];
- pnts.CopyTo(0, tpnts, 0, 3);
- CircularArc2d ca2d = GetMinCircle(tpnts[0], tpnts[1], tpnts[2], out ptlst);
- //找到点集中距离圆心的最远点为第四点
- tpnts[3] = pnts.FindByMaxKey(pnt => pnt.GetDistanceTo(ca2d.Center));
- //如果最远点属于圆结束
- while (!ca2d.IsPartOf(tpnts[3]))
- {
- //如果最远点不属于圆,按此四点计算最小包围圆
- ca2d = GetMinCircle(tpnts[0], tpnts[1], tpnts[2], tpnts[3], out ptlst);
- //将结果作为新的前三点
- if (ptlst.Count == 3)
- {
- tpnts[2] = ptlst.Last.Value;
- }
- else
- {
- //如果计算的结果只有两个点
- //选择上一步计算的圆心为基准,找出点集中圆心的最远点作为第三点
- //似乎是算法的问题?这里第三点不能任意选择,否则无法收敛
- //tpnts[2] = pnts.GetMaxBy(pnt => pnt.GetDistanceTo(ca2d.Center));
- //if (ca2d.IsPartOf(tpnts[2]))
- // return ca2d;
- //看来是理解有误,第三点应该取另两点中距离圆心较远的点
- //但按算法中描述的任选其中一点的话,还是无法收敛......
- tpnts[2] =
- tpnts.Except(ptlst)
- .FindByMaxKey(pnt => ca2d.Center.GetDistanceTo(pnt));
- }
- tpnts[0] = ptlst.First.Value;
- tpnts[1] = ptlst.First.Next.Value;
- //按此三点计算最小包围圆
- ca2d = GetMinCircle(tpnts[0], tpnts[1], tpnts[2], out ptlst);
- //找到点集中圆心的最远点为第四点
- tpnts[3] = pnts.FindByMaxKey(pnt => pnt.GetDistanceTo(ca2d.Center));
- }
- return ca2d;
- }
- public static ConvexHull2d GetConvexHull(this List<Point2d> pnts)
- {
- ConvexHull2d ch2d = new ConvexHull2d();
- ch2d.MelkmanEval(pnts);
- return ch2d;
- }
- public static ConvexHull2d GetConvexHull2(this List<Point2d> pnts)
- {
- ConvexHull2d ch2d = new ConvexHull2d();
- ch2d.JarvisEval(pnts);
- return ch2d;
- }
- public static ConvexHull2d GetConvexHull3(this List<Point2d> pnts)
- {
- ConvexHull2d ch2d = new ConvexHull2d();
- ch2d.GrahamEval(pnts);
- return ch2d;
- }
- #endregion
- }
- }
|
|