【POJ 1981】Circle and Points(已知圆上两点求圆心坐标)
【题目链接】:http://poj.org/problem?id=1981
【题意】
给你n个点(n<=300);
然后给你一个半径R;
让你在平面上找一个半径为R的圆;
这里R=1
使得这个圆覆盖的点的数目最多;
【题解】
最少会有一个点;
考虑两个点的情况;
枚举任意两个点在圆上;
考虑最极端的情况;
就是这两个点都在圆的边上;(这样圆心就尽可能地远离它们俩了,以求让这个圆覆盖更多的点);
然后求出这个时候这时的圆心的坐标;
然后看看其他的在这个圆内的点的数目就好;
圆心的话只要求一边的圆心就可以了;
不会证。
这里圆心的坐标我是用数学方法推导出来的;
可以这样
假设圆心坐标为a,b;
两个点为(x1,y1),(x2,y2)
然后根据
两式相减;
可以得到圆心所在的直线满足的方程也即得到a用b表达的形式;
然后回带入①②两个等式中任意一个;
然后解一元二次方程
得到b
从而得到a
(x1==x2的情况需要特判一下)
【Number Of WA】
0
【完整代码】
#include <iostream>
#include <cmath>
#include <algorithm>
using namespace std;
#define lson l,m,rt<<1
#define rson m+1,r,rt<<1|1
#define LL long long
#define rep1(i,a,b) for (int i = a;i <= b;i++)
#define rep2(i,a,b) for (int i = a;i >= b;i--)
#define mp make_pair
#define pb push_back
#define fi first
#define se second
#define ms(x,y) memset(x,y,sizeof x)
typedef pair<int,int> pii;
typedef pair<LL,LL> pll;
const int dx[9] = {0,1,-1,0,0,-1,-1,1,1};
const int dy[9] = {0,0,0,-1,1,-1,1,-1,1};
const double pi = acos(-1.0);
const int N = 310;
struct point
{
double x,y;
};
point t[N];
int n;
double sqr(double t)
{
return t*t;
}
double dis(point a,point b)
{
return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
point get_yx(point q,point w)
{
double R = 1.0;
double x1 = q.x,y1 = q.y;
double x2 = w.x,y2 = w.y;
double a,b,A,B;
if (fabs(x1-x2)<1e-5)
{
b = (y1+y2)/2.0;
double td = dis(q,w);
a = x1+sqrt(sqr(R)-sqr(td/2.0));
}
else
{
A=(y1-y2)/(x2-x1),B=(sqr(x2)-sqr(x1)+sqr(y2)-sqr(y1))/(2.0*(x2-x1));
double m1 = (sqr(A)+1),m2 = 2*A*B-2*y1-2*A*x1;
double m3 = sqr(x1)+sqr(B)-2*B*x1+sqr(y1)-sqr(R);
b=(-m2+sqrt(sqr(m2)-4*m1*m3))/(2*m1);
a = A*b+B;
}
point temp;
temp.x = a,temp.y = b;
return temp;
}
int main()
{
//freopen("F:\\rush.txt","r",stdin);
ios::sync_with_stdio(false);
// point dd = get_yx(point{0,1},point{1,0});
//cout << fixed << setprecision(6) <<dd.x<<' '<<dd.y<<endl;
//return 0;
while (cin>>n)
{
if (n==0) break;
rep1(i,1,n)
cin >> t[i].x>>t[i].y;
int ans = 1;
rep1(i,1,n-1)
rep1(j,i+1,n)
{
double d = dis(t[i],t[j]);
if (d>2.0) continue;
point Q = get_yx(t[i],t[j]);
int temp = 0;
rep1(k,1,n)
if (dis(t[k],Q)<1+1e-4)
temp++;
ans = max(ans,temp);
}
cout << ans << endl;
}
//printf("\n%.2lf sec \n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}