CSU 1503: 点到圆弧的距离(计算几何)

题目描述

输入一个点 P 和一条圆弧(圆周的一部分),你的任务是计算 P 到圆弧的最短距离。换句话 说,你需要在圆弧上找一个点,到 P点的距离最小。 
提示:请尽量使用精确算法。相比之下,近似算法更难通过本题的数据。 

输入

输入包含最多 10000组数据。每组数据包含 8个整数 x1, y1, x2, y2, x3, y3, xp, yp。圆弧的起点 是 A(x1,y1),经过点 B(x2,y2),结束位置是 C(x3,y3)。点 P的位置是 (xp,yp)。输入保证 A, B, C 各不相同且不会共线。上述所有点的坐标绝对值不超过 20。

输出

对于每组数据,输出测试点编号和 P 到圆弧的距离,保留三位小数。

样例输入

0 0 1 1 2 0 1 -1
3 4 0 5 -3 4 0 1 

样例输出

Case 1: 1.414
Case 2: 4.000


分两种情况:

第一种:点跟圆心的连线在那段扇形的圆弧范围内,点到圆弧的最短距离为点到圆心的距离减去半径然后取绝对值;

第二种:点跟圆心的连线不在那段扇形的圆弧范围内,点到圆弧的最短的距离为到这段圆弧的两个端点的最小值。

 

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 const double eps  = 1e-8;
 4 const double PI = acos(-1.0);
 5 int casee = 0;
 6 int dcmp(double x) {
 7     if(fabs(x) < eps) return 0;
 8     else return x < 0 ? -1 : 1;
 9 }
10 struct Vector {
11     double x, y;
12     Vector(double x_ = 0, double y_ = 0) : x(x_), y(y_) {}
13 
14     Vector operator+ (const Vector &a) const {
15         return Vector(x + a.x, y + a.y);
16     }
17     Vector operator- (const Vector &a) const {
18         return Vector(x - a.x, y - a.y);
19     }
20     Vector operator* (double p) {
21         return Vector(x * p, y * p);
22     }
23     friend Vector operator* (double p, const Vector &a) {
24         return Vector(a.x * p, a.y * p);
25     }
26     Vector operator/ (double p) {
27         return Vector(x / p, y / p);
28     }
29     bool operator== (const Vector &a) const {
30         return dcmp(x - a.x) == 0 && dcmp(y - a.y) == 0;
31     }
32 
33     friend double dmul(const Vector &a, const Vector &b) {
34         return a.x * b.x + a.y * b.y;
35     }
36     friend double cmul(const Vector &a, const Vector &b) {
37         return a.x * b.y - a.y * b.x;
38     }
39     friend double angle(const Vector &a, const Vector &b) {
40         return acos(dmul(a, b) / a.length() / b.length());
41     }
42     friend Vector rotate(const Vector &a, double rad) {
43         return Vector(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad));
44     }
45     friend Vector normal(const Vector &a) {
46         double l = a.length();
47         return Vector(a.x / l, a.y / l);
48     }
49 
50     double length() const {
51         return sqrt(dmul(*this, *this));
52     }
53 };
54 typedef Vector Point;
55 double dis(Point A, Point B) {
56     return (A-B).length();
57 }
58 Point pp1, pp2, pp3, pp, o;
59 
60 Point get_cir(Point a,Point b,Point c) {//已知圆上3点,求圆的外心
61     double a1=b.x-a.x,b1=b.y-a.y,c1=(a1*a1+b1*b1)/2;
62     double a2=c.x-a.x,b2=c.y-a.y,c2=(a2*a2+b2*b2)/2;
63     double d=a1*b2-a2*b1;
64     return Point(a.x+(c1*b2-c2*b1)/d,a.y+(a1*c2-a2*c1)/d);
65 }
66 double AngleToX(Vector A, Vector B) {//向量与x轴正方向的夹角,范围[0,2*pi) B为(1,0)
67     double res = angle(A, B);
68     if(dcmp(cmul(A, B)) < 0) res = PI * 2 - res;
69     return res;
70 }
71 void gao (Point A,Point B,Point C,Point P,Point O)
72 {
73     Vector X(1, 0);
74     Vector OA = A - O, OP = P - O, OC = C - O, OB = B - O;
75     double fa = AngleToX(OA, X);
76     double fb = AngleToX(OB, X);
77     double fc = AngleToX(OC, X);
78     double fp = AngleToX(OP, X);
79     double ans = min(dis(P, A), dis(P, C));
80     double r = dis(A, O);
81     if(fa>fc) swap(fa,fc);
82     if(dcmp(fa-fb)<=0 && dcmp(fb - fc) <= 0 && dcmp(fa - fp) <= 0 && dcmp(fp - fc) <= 0) {
83         double r = dis(A, O);
84         ans = min(ans, fabs(dis(O, P) - r));
85     }
86     if(!(dcmp(fa - fb) <= 0 &&dcmp(fb-fc)<= 0)&&!(dcmp(fa-fp)<=0&&dcmp(fp-fc)<=0)) {
87         ans = min(ans, fabs(dis(O, P)- r));
88     }
89 
90     printf("Case %d: %.3f\n", ++casee, ans);
91 }
92 int main() {
93     while(~scanf("%lf%lf%lf%lf%lf%lf%lf%lf", &pp1.x,&pp1.y,&pp2.x,&pp2.y,&pp3.x,&pp3.y,&pp.x,&pp.y)){
94         o=get_cir(pp1,pp2,pp3);
95         gao(pp1,pp2,pp3,pp,o);
96     }
97     return 0;
98 }

解法2

 1 #include <cstdio>
 2 #include <cstring>
 3 #include <algorithm>
 4 #include <cmath>
 5 #include <map>
 6 #define eps (1e-4)
 7 #define N 205
 8 #define dd double
 9 #define sqr(x) ((x)*(x))
10 const double pi = acos(-1);
11 using namespace std;
12 struct Point
13 {
14     double x,y;
15 };
16 double cross(Point a,Point b,Point c) ///叉积
17 {
18     return (a.x-c.x)*(b.y-c.y)-(b.x-c.x)*(a.y-c.y);
19 }
20 double dis(Point a,Point b) ///距离
21 {
22     return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
23 }
24 Point waixin(Point a,Point b,Point c) ///外接圆圆心坐标
25 {
26     Point p;
27     double a1 = b.x - a.x, b1 = b.y - a.y, c1 = (a1*a1 + b1*b1)/2;
28     double a2 = c.x - a.x, b2 = c.y - a.y, c2 = (a2*a2 + b2*b2)/2;
29     double d = a1*b2 - a2*b1;
30     p.x = a.x + (c1*b2 - c2*b1)/d, p.y=a.y + (a1*c2 -a2*c1)/d;
31     return p;
32 }
33 bool judge(Point a,Point b,Point c,Point p){ ///此模板判断点 p 是否在两条射线 ab 和 ac 之间(但是p点不在 ab或者 ac上)
34     if(cross(b,c,a)>0&&cross(p,a,c)<0&&cross(p,a,b)>0){ ///b 在 c 的逆时针方向
35         return true;
36     }
37     if(cross(b,c,a)<0&&cross(p,a,c)>0&&cross(p,a,b)<0){ ///b 在 c 的顺时针方向
38         return true;
39     }
40     return false;
41 }
42 bool judge1(Point a,Point b,Point c,Point p){ ///此模板判断点 p 是否在两条射线 ab 和 ac 之间(p点可以在 ab或者 ac上)
43     if(cross(b,c,a)>0&&cross(p,a,c)<=0&&cross(p,a,b)>=0){ ///b 在 c 的逆时针方向
44         return true;
45     }
46     if(cross(b,c,a)<0&&cross(p,a,c)>=0&&cross(p,a,b)<=0){ ///b 在 c 的顺时针方向
47         return true;
48     }
49     return false;
50 }
51 int main()
52 {
53     Point a,b,c,p;
54     int t = 1;
55     freopen("a.in","r",stdin);
56     freopen("a.txt","w",stdout);
57     while(scanf("%lf%lf%lf%lf%lf%lf%lf%lf",&a.x,&a.y,&b.x,&b.y,&c.x,&c.y,&p.x,&p.y)!=EOF)
58     {
59 
60         Point circle = waixin(a,b,c);
61         double r = dis(circle,a);
62         double ans = min(dis(p,a),dis(p,c));
63         double op = dis(circle,p);
64         if(fabs(2*r-dis(a,c))<eps){ ///平角特殊处理
65             if(cross(p,c,circle)*cross(b,c,circle)>=0){
66                 ans = min(ans,fabs(op-r));
67             }
68             printf("Case %d: %.3lf\n",t++,ans);
69             continue;
70         }
71         if(judge(circle,a,c,b))  ///劣弧
72         {
73             if(judge1(circle,a,c,p)) ans = min(ans,fabs(op-r));
74         }
75         else
76         {
77             if(!judge(circle,a,c,p)) ans = min(ans,fabs(op-r));
78         }
79         printf("Case %d: %.3lf\n",t++,ans);
80     }
81 }

标程:

 1 // Rujia Liu
 2 #include<cmath>
 3 #include<cstdio>
 4 #include<iostream>
 5 
 6 using namespace std;
 7 
 8 const double PI = acos(-1.0);
 9 const double TWO_PI = 2 * PI;
10 const double eps = 1e-6;
11 
12 inline double NormalizeAngle(double rad, double center = PI) {
13   return rad - TWO_PI * floor((rad + PI - center) / TWO_PI);
14 }
15 
16 inline int dcmp(double x) {
17   if(fabs(x) < eps) return 0; else return x < 0 ? -1 : 1;
18 }
19 
20 struct Point {
21   double x, y;
22   Point(double x=0, double y=0):x(x),y(y) { }
23 };
24 
25 typedef Point Vector;
26 
27 inline Vector operator + (Vector A, Vector B) { return Vector(A.x+B.x, A.y+B.y); }
28 inline Vector operator - (Point A, Point B) { return Vector(A.x-B.x, A.y-B.y); }
29 
30 inline double Dot(Vector A, Vector B) { return A.x*B.x + A.y*B.y; }
31 inline double Cross(Vector A, Vector B) { return A.x*B.y - A.y*B.x; }
32 inline double Length(Vector A) { return sqrt(Dot(A, A)); }
33 
34 // 外接圆圆心。假定三点不共线
35 Point get_circumscribed_center(Point p1, Point p2, Point p3) {
36   double bx = p2.x - p1.x;
37   double by = p2.y - p1.y;
38   double cx = p3.x - p1.x;
39   double cy = p3.y - p1.y;
40   double d = 2 * (bx * cy - by * cx);
41   Point p;
42   p.x = (cy * (bx * bx + by * by) - by * (cx * cx + cy * cy)) / d + p1.x;
43   p.y = (bx * (cx * cx + cy * cy) - cx * (bx * bx + by * by)) / d + p1.y;
44   return p;
45 }
46 
47 double DistanceToArc(Point a, Point start, Point mid, Point end) {   ///点到圆弧的距离
48   Point p = get_circumscribed_center(start, mid, end);
49   bool CCW = dcmp(Cross(mid - start, end - start)) > 0;
50   double ang_start = atan2(start.y-p.y, start.x-p.x);
51   double ang_end = atan2(end.y-p.y, end.x-p.x);
52   double r = Length(p - start);
53   double ang = atan2(a.y-p.y, a.x-p.x);
54   bool inside;
55   if(CCW) {
56     inside = NormalizeAngle(ang - ang_start) < NormalizeAngle(ang_end - ang_start);
57   } else {
58     inside = NormalizeAngle(ang - ang_end) < NormalizeAngle(ang_start - ang_end);
59   }
60   if(inside) {
61     return fabs(r - Length(p - a));
62   }
63   return min(Length(a - start), Length(a - end));
64 }
65 
66 int main() {
67   int kase = 0;
68   double x1, y1, x2, y2, x3, y3, xp, yp;
69   while(cin >> x1 >> y1 >> x2 >> y2 >> x3 >> y3 >> xp >> yp) {
70     double ans = DistanceToArc(Point(xp,yp), Point(x1,y1), Point(x2,y2), Point(x3,y3));
71     printf("Case %d: %.3lf\n", ++kase, ans);
72   }
73   return 0;
74 }

 

转载于:https://www.cnblogs.com/agenthtb/p/7635259.html

  • 0
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值