• 最小覆盖圆算法


    最小圆覆盖,很经典的问题。题目大概是,平面上n个点,求一个半径最小的圆,能够覆盖所有的点。

     

    算法有点难懂,于是讲讲我的理解。

    如果要求一个最小覆盖圆,这个圆至少要由三个点确定。有一种算法就是任意取三个点作圆,然后判断距离圆心最远的点是否在圆内,若在,则完成;若不在则用最远点更新这个圆。

    这里介绍的算法是,先任意选取两个点,以这两个点的连线为直径作圆。再以此判断剩余的点,看它们是否都在圆内(或圆上),如果都在,说明这个圆已经找到。如果没有都在:假设我们用的最开始的两个点为p[1],p[2],并且找到的第一个不在圆内(或圆上)的点为p[i],于是我们用这个点p[i]去寻找覆盖p[1]到p[i-1]的最小覆盖圆。

    那么,过确定点p[i]的从p[1]到p[i-1]的最小覆盖圆应该如何求呢?

    我们先用p[1]和p[i]做圆,再从2到i-1判断是否有点不在这个圆上,如果都在,则说明已经找到覆盖1到i-1的圆。如果没有都在:假设我们找到第一个不在这个圆上的点为p[j],于是我们用两个已知点p[j]与p[i]去找覆盖1到j-1的最小覆盖圆。

    而对于两个已知点p[j]与p[i]求最小覆盖圆,只要从1到j-1中,第k个点求过p[k],p[j],p[i]三个点的圆,再判断k+1到j-1是否都在圆上,若都在,说明找到圆;若有不在的,则再用新的点p[k]更新圆即可。

    于是,这个问题就被转化为若干个子问题来求解了。

    由于三个点确定一个圆,我们的过程大致上做的是从没有确定点,到有一个确定点,再到有两个确定点,再到有三个确定点来求圆的工作。

    关于正确性的证明以及复杂度的计算这里就不介绍了,可以去看完整的算法介绍:http://wenku.baidu.com/view/162699d63186bceb19e8bbe6.html

    恩。关于细节方面。

    a.通过三个点如何求圆?

       先求叉积。

       若叉积为0,即三个点在同一直线,那么找到距离最远的一对点,以它们的连线为直径做圆即可;

       若叉积不为0,即三个点不共线,那么就是第二个问题,如何求三角形的外接圆?

    b.如何求三角形外接圆?

       假设三个点(x1,y1),(x2,y2),(x3,y3);

       设过(x1,y1),(x2,y2)的直线l1方程为Ax+By=C,它的中点为(midx,midy)=((x1+x2)/2,(y1+y2)/2),l1中垂线方程为A1x+B1y=C1;则它的中垂线方程中A1=-B=x2-x1,B1=A=y2-y1,C1=-B*midx+A*midy=((x2^2-x1^2)+(y2^2-y1^2))/2;

       同理可以知道过(x1,y1),(x3,y3)的直线的中垂线的方程。

       于是这两条中垂线的交点就是圆心。

    c.如何求两条直线交点?

       设两条直线为A1x+B1y=C1和A2x+B2y=C2。

       设一个变量det=A1*B2-A2*B1;

       如果det=0,说明两直线平行;若不等于0,则求交点:x=(B2*C1 -B1*C2)/det,y=(A1*C2-A2*C1)/det;

    d.于是木有了。。

     1 #include<stdio.h>
     2 #include<math.h>
     3 struct   TPoint  
     4 {  
     5          double x,y;  
     6 };
     7 TPoint a[1005],d;
     8 double r;
     9 
    10 double   distance(TPoint   p1,   TPoint   p2)   //两点间距离
    11 {  
    12     return (sqrt((p1.x-p2.x)*(p1.x -p2.x)+(p1.y-p2.y)*(p1.y-p2.y)));      
    13 }
    14 double multiply(TPoint   p1,   TPoint   p2,   TPoint   p0)  
    15 {  
    16     return   ((p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y));          
    17 }  
    18 void MiniDiscWith2Point(TPoint   p,TPoint   q,int   n)
    19 {
    20     d.x=(p.x+q.x)/2.0;
    21     d.y=(p.y+q.y)/2.0;
    22     r=distance(p,q)/2;
    23     int k;
    24     double c1,c2,t1,t2,t3;
    25     for(k=1;k<=n;k++)
    26     {
    27         if(distance(d,a[k])<=r)continue;
    28         if(multiply(p,q,a[k])!=0.0)
    29         {
    30             c1=(p.x*p.x+p.y*p.y-q.x*q.x-q.y*q.y)/2.0;
    31             c2=(p.x*p.x+p.y*p.y-a[k].x*a[k].x-a[k].y*a[k].y)/2.0;
    32 
    33             d.x=(c1*(p.y-a[k].y)-c2*(p.y-q.y))/((p.x-q.x)*(p.y-a[k].y)-(p.x-a[k].x)*(p.y-q.y));
    34             d.y=(c1*(p.x-a[k].x)-c2*(p.x-q.x))/((p.y-q.y)*(p.x-a[k].x)-(p.y-a[k].y)*(p.x-q.x));
    35             r=distance(d,a[k]);
    36         }
    37         else
    38         {
    39             t1=distance(p,q);
    40             t2=distance(q,a[k]);
    41             t3=distance(p,a[k]);
    42             if(t1>=t2&&t1>=t3)
    43             {d.x=(p.x+q.x)/2.0;d.y=(p.y+q.y)/2.0;r=distance(p,q)/2.0;}
    44             else if(t2>=t1&&t2>=t3)
    45             {d.x=(a[k].x+q.x)/2.0;d.y=(a[k].y+q.y)/2.0;r=distance(a[k],q)/2.0;}
    46             else
    47             {d.x=(a[k].x+p.x)/2.0;d.y=(a[k].y+p.y)/2.0;r=distance(a[k],p)/2.0;}
    48         }
    49     }
    50 }
    51 
    52 void MiniDiscWithPoint(TPoint   pi,int   n)
    53 {
    54     d.x=(pi.x+a[1].x)/2.0;
    55     d.y=(pi.y+a[1].y)/2.0;
    56     r=distance(pi,a[1])/2.0;
    57     int j;
    58     for(j=2;j<=n;j++)
    59     {
    60     if(distance(d,a[j])<=r)continue;
    61         else
    62         {
    63             MiniDiscWith2Point(pi,a[j],j-1);
    64         }
    65     }
    66 }
    67 int main()
    68 {
    69     int i,n;
    70     while(scanf("%d",&n)&&n)
    71     {
    72         for(i=1;i<=n;i++)
    73         {
    74             scanf("%lf %lf",&a[i].x,&a[i].y);
    75         }
    76         if(n==1)
    77         { printf("%.2lf %.2lf 0.00
    ",a[1].x,a[1].y);continue;}
    78         r=distance(a[1],a[2])/2.0;
    79         d.x=(a[1].x+a[2].x)/2.0;
    80         d.y=(a[1].y+a[2].y)/2.0;
    81         for(i=3;i<=n;i++)
    82         {
    83             if(distance(d,a[i])<=r)continue;
    84             else
    85             MiniDiscWithPoint(a[i],i-1);
    86         }
    87         printf("%.2lf %.2lf %.2lf
    ",d.x,d.y,r);
    88     }
    89     return 0;
    90 }
    View Code
  • 相关阅读:
    远程桌面连接win10问题解决
    为什么n各节点的的二叉链表中有n+1个空链域
    西门子Step7找不到有效授权的解决方法
    表达式树获取函数命名
    逆波兰表达式
    双向循环链表实践
    快速找到未知长度单链表的中间节点
    java的ArrayList(线性表)和LinkedList(双向链表)的深入学习
    23种设计模式中的访问者模式
    23种设计模式中的原型模式
  • 原文地址:https://www.cnblogs.com/lxm940130740/p/3900871.html
Copyright © 2020-2023  润新知