zoukankan      html  css  js  c++  java
  • wenbao与最小圆覆盖

    一下来自大神博客:http://blog.sina.com.cn/s/blog_6e63f59e010120dl.html

    这里介绍的算法是,先任意选取两个点,以这两个点的连线为直径作圆。再以此判断剩余的点,看它们是否都在圆内(或圆上),如果都在,说明这个圆已经找到。如果没有都在:假设我们用的最开始的两个点为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]更新圆即可。

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

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

    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.于是木有了。。

    http://oj.jxust.edu.cn/problem.php?id=1328

     大牛代码:

     1 #include<stdio.h>
     2 #include<math.h>
     3 struct TPoint{  
     4     double x,y;  
     5 };
     6 TPoint a[1005],d;
     7 double r;
     8 
     9 double   distance(TPoint   p1,   TPoint   p2){  
    10     return (sqrt((p1.x-p2.x)*(p1.x -p2.x)+(p1.y-p2.y)*(p1.y-p2.y)));      
    11 }
    12 double multiply(TPoint   p1,   TPoint   p2,   TPoint   p0)  {  
    13     return   ((p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y));          
    14 }  
    15 void MiniDiscWith2Point(TPoint p,TPoint q,int n){
    16     d.x=(p.x+q.x)/2.0;
    17     d.y=(p.y+q.y)/2.0;
    18     r=distance(p,q)/2;
    19     int k;
    20     double c1,c2,t1,t2,t3;
    21     for(k=1;k<=n;k++){
    22         if(distance(d,a[k])<=r)continue;
    23         if(multiply(p,q,a[k])!=0.0){
    24             c1=(p.x*p.x+p.y*p.y-q.x*q.x-q.y*q.y)/2.0;
    25             c2=(p.x*p.x+p.y*p.y-a[k].x*a[k].x-a[k].y*a[k].y)/2.0;
    26             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));
    27             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));
    28             r=distance(d,a[k]);
    29         }
    30         else{
    31             t1=distance(p,q);
    32             t2=distance(q,a[k]);
    33             t3=distance(p,a[k]);
    34             if(t1>=t2&&t1>=t3)
    35             {d.x=(p.x+q.x)/2.0;d.y=(p.y+q.y)/2.0;r=distance(p,q)/2.0;}
    36             else if(t2>=t1&&t2>=t3)
    37             {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;}
    38             else
    39             {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;}
    40         }
    41     }
    42 }
    43 
    44 void MiniDiscWithPoint(TPoint   pi,int   n){
    45     d.x=(pi.x+a[1].x)/2.0;
    46     d.y=(pi.y+a[1].y)/2.0;
    47     r=distance(pi,a[1])/2.0;
    48     int j;
    49     for(j=2;j<=n;j++){
    50         if(distance(d,a[j])<=r)continue;
    51         else{
    52             MiniDiscWith2Point(pi,a[j],j-1);
    53         }
    54     }
    55 }
    56 int main(){
    57     int i,n, t;
    58     scanf("%d ", &t);
    59     while(t--){
    60         scanf("%d",&n);
    61         for(i=1;i<=n;i++){
    62             scanf("%lf %lf",&a[i].x,&a[i].y);
    63         }
    64         if(n==1)
    65         { printf("%.1lf %.1lf
    ",a[1].x,a[1].y);continue;}
    66         r=distance(a[1],a[2])/2.0;
    67         d.x=(a[1].x+a[2].x)/2.0;
    68         d.y=(a[1].y+a[2].y)/2.0;
    69         for(i=3;i<=n;i++){
    70             if(distance(d,a[i])<=r)continue;
    71             else
    72                 MiniDiscWithPoint(a[i],i-1);
    73         }
    74         printf("%.1lf %.1lf
    ",d.x,d.y);
    75     }
    76     return 0;
    77 }

     垃圾代码:

     1 #include <iostream>
     2 #include <stdio.h>
     3 #include <cmath>
     4 using namespace std;
     5 double r;
     6 struct Node {
     7     double x, y;
     8 } T[1005], t;
     9 double dis(Node a, Node b){
    10     return sqrt((a.x - b.x)*(a.x - b.x) + (a.y - b.y)*(a.y - b.y));
    11 }
    12 double xmup(Node a, Node b, Node c){
    13     return (a.x - b.x)*(a.y - c.y) - (a.y - b.y)*(a.x - c.x);
    14 }
    15 void FindMin2(Node a, Node b, int x){
    16     t.x = (a.x + b.x) /2.0;
    17     t.y = (a.y + b.y) /2.0;
    18     r = dis(a, b) /2.0;
    19     for(int i = 1; i < x; i++){
    20         if(dis(t, T[i]) <= r) continue;
    21         if(xmup(a, b, T[i]) == 0.0){
    22             t.x = (a.x + b.x + T[i].x) /3.0;
    23             t.y = (a.y + b.y + T[i].y) /3.0;
    24             r = dis(t, T[i]);
    25         }else{
    26             double a1 = 2.0*(a.x - b.x), a2 = 2.0*(a.x - T[i].x);
    27             double b1 = 2.0*(a.y - b.y), b2 = 2.0*(a.y - T[i].y);
    28             double c1 = a.x*a.x + a.y*a.y - b.x*b.x - b.y*b.y;
    29             double c2 = a.x*a.x + a.y*a.y - T[i].x*T[i].x - T[i].y*T[i].y;
    30             t.x = (c1*b2 - c2*b1) / (a1*b2 - b1*a2);
    31             t.y = (a1*c2 - a2*c1) / (a1*b2 - b1*a2);
    32             r = dis(t, T[i]);
    33         }
    34     }
    35 }
    36 void FindMin(Node a, int x){
    37     t.x = (a.x + T[1].x) /2.0;
    38     t.y = (a.y + T[1].y) /2.0;
    39     r = dis(a, T[1]) /2.0;
    40     for(int i = 2; i < x; i++){
    41         if(dis(t, T[i]) <= r) continue;
    42         FindMin2(a, T[i], i);
    43     }
    44 }
    45 int main(){
    46     int n, m;
    47     scanf("%d", &n);
    48     while(n--){
    49         scanf("%d", &m);
    50         for(int i = 1; i <= m; i++){
    51             scanf("%lf %lf", &T[i].x, &T[i].y);
    52         }
    53         if(m == 1){
    54             printf("%.1lf %.1lf
    ", T[1].x, T[1].y);
    55             continue;
    56         }
    57         t.x = (T[1].x + T[2].x) /2.0;
    58         t.y = (T[1].y + T[2].y) /2.0;
    59         r = dis(T[1], T[2]) /2.0;
    60         for(int i = 3; i <= m; i++){
    61             if(dis(t, T[i]) <= r) continue;
    62             FindMin(T[i], i);
    63         }
    64         printf("%.1lf %.1lf
    ", t.x, t.y);
    65     }
    66     return 0;
    67 }

    http://acm.hdu.edu.cn/showproblem.php?pid=3007

     1 #include <cmath>
     2 #include <stdio.h>
     3 #include <iostream>
     4 using namespace std;
     5 double r;
     6 struct Node{
     7     double x, y;
     8 }T[505], t;
     9 double dis(Node a, Node b){
    10     return sqrt((a.x - b.x)*(a.x - b.x) + (a.y - b.y)*(a.y-b.y));
    11 }
    12 double xmup(Node a, Node b, Node c){
    13     return (a.x - b.x)*(c.y - b.y) - (a.y - b.y)*(c.x - b.x);
    14 }
    15 void FindMin2(Node a, Node b, int x){
    16     t.x = (a.x + b.x) /2.0;
    17     t.y = (a.y + b.y) /2.0;
    18     r = dis(a, b) /2.0;
    19     for(int i = 1; i < x; i++){
    20         if(dis(t, T[i]) <= r) continue;
    21         if(xmup(a, b, T[i]) == 0.0) {
    22             t.x = (a.x + b.x + T[i].x) /3.0;
    23             t.y = (a.y + b.y + T[i].y) /3.0;
    24             r = dis(t, T[i]);
    25         }else{
    26             double a1 = 2.0*(a.x - b.x), a2 = 2.0*(a.x - T[i].x);
    27             double b1 = 2.0*(a.y - b.y), b2 = 2.0*(a.y - T[i].y);
    28             double c1 = a.x*a.x + a.y*a.y - b.x*b.x - b.y*b.y;
    29             double c2 = a.x*a.x + a.y*a.y - T[i].x*T[i].x - T[i].y*T[i].y;
    30             t.x = (c1*b2 - c2*b1) / (a1*b2 - b1*a2);
    31             t.y = (a1*c2 - a2*c1) / (a1*b2 - b1*a2);
    32             r = dis(t, T[i]);
    33         }
    34     }
    35 }
    36 void FindMin(Node a, int x){
    37     t.x = (T[1].x + a.x) /2.0;
    38     t.y = (T[1].y + a.y) /2.0;
    39     r = dis(T[1], a) /2.0;
    40     for(int i = 2; i < x; i++){
    41         if(dis(t, T[i]) <= r) continue;
    42         else FindMin2(a, T[i], i);
    43     }
    44 }
    45 int main(){
    46     int n;
    47     while(scanf("%d", &n)){
    48         if(n == 0) break;
    49         for(int i = 1; i <= n; i++){
    50             scanf("%lf %lf", &T[i].x, &T[i].y);
    51         }
    52         if(n == 1){
    53             printf("wenbaozuishuai");
    54             continue;
    55         }
    56         t.x = (T[1].x + T[2].x) /2.0;
    57         t.y = (T[1].y + T[2].y) /2.0;
    58         r = dis(T[1], T[2]) /2.0;
    59         for(int i = 3; i <= n; i++){
    60             if(dis(t, T[i]) <= r) continue;
    61             FindMin(T[i], i);
    62         }
    63         printf("%.2lf %.2lf %.2lf
    ", t.x, t.y, r);
    64     }
    65     return 0;
    66 }

    只有不断学习才能进步!

  • 相关阅读:
    centos 用户管理
    rsync 实验
    文件共享和传输
    PAT 1109 Group Photo
    PAT 1108 Finding Average
    PAT 1107 Social Clusters
    PAT 1106 Lowest Price in Supply Chain
    PAT 1105 Spiral Matrix
    PAT 1104 Sum of Number Segments
    PAT 1103 Integer Factorization
  • 原文地址:https://www.cnblogs.com/wenbao/p/6015541.html
Copyright © 2011-2022 走看看