Go: package main

import (
	"fmt"
	"math"
	"math/cmplx"
)

func main() {
	for {
		var n int
		fmt.Scan(&n)

		if n == 0 {
			return
		}

		arr := make([]complex128, n)

		for i := 0; i < n; i++ {
			var x, y float64

			fmt.Scanf("%f%f", &x, &y)

			arr[i] = complex(x, y)
		}

		max := 1
		for i := 0; i < n; i++ {
			for j := i + 1; j < n; j++ {
				sub := arr[j] - arr[i]

				if cmplx.Abs(sub) >= 2 {
					continue
				}
				sum := arr[j] + arr[i]

				v := 1
				{
					halfd := cmplx.Abs(sub / 2)

					center := sum/2 + sub*cmplx.Rect(math.Sqrt(1-halfd*halfd)/cmplx.Abs(sub), float64(v)*math.Pi/2)

					counter := 0
					for k := 0; k < n; k++ {
						if i == k || j == k || cmplx.Abs(arr[k]-center) <= 1 {
							counter++
						}
					}
					if max < counter {
						max = counter
					}
				}

				v = -1
				{
					halfd := cmplx.Abs(sub / 2)

					center := sum/2 + sub*cmplx.Rect(math.Sqrt(1-halfd*halfd)/cmplx.Abs(sub), float64(v)*math.Pi/2)

					counter := 0
					for k := 0; k < n; k++ {
						if i == k || j == k || cmplx.Abs(arr[k]-center) <= 1 {
							counter++
						}
					}
					if max < counter {
						max = counter
					}
				}

			}
		}

		fmt.Println(max)
	}
}

C++: #include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef pair<int,int> P;
typedef pair<int,P> P1;
typedef pair<P,P> P2;
#define pu push
#define pb push_back
#define mp make_pair
#define INF 1000000000
#define fi first
#define sc second
#define rep(i,x) for(int i=0;i<x;i++)
#define SORT(x) sort(x.begin(),x.end())
#define ERASE(x) x.erase(unique(x.begin(),x.end()),x.end())
#define POSL(x,v) (lower_bound(x.begin(),x.end(),v)-x.begin())
#define POSU(x,v) (upper_bound(x.begin(),x.end(),v)-x.begin())
typedef complex<double> pt;
typedef pair<pt,pt> L;
typedef vector<P> poly;
const double EPS = 1e-15;
#define x real()
#define y imag()
 
double dot(pt a,pt b){
	return (conj(a)*b).x;
}
double cross(pt a,pt b){
	return (conj(a)*b).y;
}
int ccw(pt a,pt b,pt c){
	b -= a; c -= a;
	if(cross(b,c) > EPS) return 1; // counter clockwise
	if(cross(b,c) < -EPS) return -1; // clockwise
	if(dot(b,c) < -EPS) return 2; //c-a-b
	if(norm(b) < norm(c)) return -2; //a-b-c
	return 0; //a-c-b
}
bool cmp(const pt& a,const pt& b){
	if(-EPS < a.x-b.x && a.x-b.x < EPS) return a.y < b.y;
	else return a.x < b.x;
}
vector<pt>convex_hull(vector<pt>ps)
{
	sort(ps.begin(),ps.end(),cmp);
	int k=0,n = ps.size();
	vector<pt>qs(n*2);
	for(int i=0;i<n;i++)
	{
		while(k>1 && ccw(qs[k-2],qs[k-1],ps[i]) == -1) k--;
		qs[k++] = ps[i];
	}
	for(int i=n-2,t=k;i>=0;i--)
	{
		while(k>t && ccw(qs[k-2],qs[k-1],ps[i]) == -1) k--;
		qs[k++] = ps[i];
	}
	qs.resize(k-1);
	return qs;
}
pt intersect(pt a,pt b,pt c,pt d){
	double dev = (b.y-a.y)*(d.x-c.x)-(b.x-a.x)*(d.y-c.y);
	
	double d1, d2;
 
	d1 = (c.y*d.x-c.x*d.y);
	d2 = (a.y*b.x-a.x*b.y);
 
	pt ap1;
	ap1.x = d1*(b.x-a.x) - d2*(d.x-c.x);
	ap1.x /= dev;
	ap1.y = d1*(b.y-a.y) - d2*(d.y-c.y);
	ap1.y /= dev;
 
	return ap1;
}
bool par(pt a,pt b,pt c,pt d){
	b.x -= a.x; b.y -= a.y;
	d.x -= c.x; d.y -= c.y;
	if(abs(b.y*d.x - b.x*d.y) < EPS){
		return 1;
	}
	return 0;
}
double PI = 3.1415926535897932;
int main(){
	while(1){
		int n; scanf("%d",&n); if(!n) return 0; 
		pt p[305];
		for(int i=1;i<=n;i++){
			double a,b; scanf("%lf%lf",&a,&b);
			p[i] = pt(a,b);
		}
		int ans = 1;
		for(int i=1;i<=n;i++){
			for(int j=i+1;j<=n;j++){
				pt q = (p[i]+p[j])/2.0;
				pt r = (p[i]-p[j]); if(abs(r) >= 2.0) continue;
				pt p1 = q + r / abs(r) * polar(sqrt(1.0-abs(r)*abs(r)/4.0),PI/2);
				pt p2 = q - r / abs(r) * polar(sqrt(1.0-abs(r)*abs(r)/4.0),PI/2);
				int c=0,d=0;
				for(int k=1;k<=n;k++){
					if(k==i || k==j || abs(p1-p[k]) <= 1.0+EPS) c++;
					if(k==i || k==j || abs(p2-p[k]) <= 1.0+EPS) d++;
				}
				ans = max(ans,max(c,d));
			}
		}
		printf("%d\n",ans);
	}
}

C++: #include <bits/stdc++.h>
using namespace std;
#define int long long   // <-----!!!!!!!!!!!!!!!!!!!

#define rep(i,n) for (int i=0;i<(n);i++)
#define rep2(i,a,b) for (int i=(a);i<(b);i++)
#define rrep(i,n) for (int i=(n)-1;i>=0;i--)
#define rrep2(i,a,b) for (int i=(a)-1;i>=b;i--)
#define all(a) (a).begin(),(a).end()
#define rall(a) (a).rbegin(),(a).rend()
#define printV(_v) for(auto _x:_v){cout<<_x<<" ";}cout<<endl
#define printVS(vs) for(auto x : vs){cout << x << endl;}
#define printVV(_vv) for(auto _v:_vv){for(auto _x:_v){cout<<_x<<" ";}cout<<endl;}
#define printP(p) cout << p.first << " " << p.second << endl
#define printVP(vp) for(auto p : vp) printP(p);

typedef long long ll;
typedef pair<int, int> Pii;
typedef tuple<int, int, int> TUPLE;
typedef vector<int> vi;
typedef vector<vi> vvi;
typedef vector<vvi> vvvi;
typedef vector<Pii> vp;
const int inf = 1e9;
const int mod = 1e9 + 7;

// #define double long double

typedef complex<double> P;
typedef vector<P> G;
#define here(g, i) g[i]
#define next(g, i) g[(i + 1) % g.size()]
#define prev(g, i) g[(i - 1 + g.size()) % g.size()]
const double EPS = 1e-6;
const double INF = 1e12;
const double PI = acos(-1);

P readP() {
    double x, y;
    cin >> x >> y;
    return P(x, y);
}

// rot p around q by theta (counter-clockwise)
P rotP(P p, P q, double theta) {
    p -= q;
    double x = p.real(), y = p.imag();
    p = P(x * cos(theta) - y * sin(theta), x * sin(theta) + y * cos(theta));
    p += q;
    return p;
}

P center(P p, P q) {
//     P v = q - p;
//     double theta = asin(abs(v) / 2);
//     v *= (1 / abs(v));
//     P u = rotP(v, p, theta);
//     return p + u;

    P v = q - p;
    P h = v * P(0., 1.) / abs(v);
    return p + (v / 2.) + h;
}

signed main() {
    std::ios::sync_with_stdio(false);
    std::cin.tie(0);

    int n;
    while (cin >> n, n) {
        vector<P> p(n);
        rep(i, n) p[i] = readP();

        int ans = 1;
        rep(i, n) {
            rep(j, n) {
                if (i == j) continue;
                if (norm(p[i] - p[j]) > 4. + EPS) continue;
                P c = center(p[i], p[j]);
                int cnt = 2;
                rep(k, n) {
                    if (k == i || k == j) continue;
                    if (norm(p[k] - c) < 1. + EPS) {
                        cnt++;
                    }
                }
                ans = max(ans, cnt);
            }
        }

        cout << ans << endl;
    }
}