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 "stdio.h"
#include "math.h"

#include <vector>
#include <map>
#include <set>
#include <stack>
#include <queue>
#include <cstring>
#include <string>
#include <sstream>
#include <algorithm>
#include <iomanip>
#include <iostream>

#define VARIABLE(x) cerr << #x << "=" << x << endl
#define rep(i,n) for(int i=0;i<(int)(n);i++)
#define REP(i,m,n) for (int i=m;i<(int)(n);i++)

const int INF = 10000000;
int dx[4]={1, 0, -1, 0}, dy[4]={0, 1, 0, -1};
using namespace std;
typedef long long  ll;
typedef pair<double, double> P;

/** Problem1132 : Circle and Point **/
int main()
{
	int N, ans;
	double X, Y;
	vector<P> points;
	
	while (cin>>N, N) {
		ans=0;
		
		points.clear();
		rep(i, N) {
			cin >> X >> Y;
			points.push_back(P(X, Y));
		}

		for (double x=0.0; x<=10.0; x+=0.01) {
			for (double y=0.0; y<=10.0; y+=0.01) {
				int num=0;
				vector<P>::iterator it = points.begin();
				for (; it!=points.end(); it++) {
					double d = (it->first-x)*(it->first-x) + (it->second-y)*(it->second-y);
					if (d<1.001) {
						num++;
					}
				}
				if (ans<num)
					ans = num;
			}
		}
		cout << ans << endl;
	}
	return 0;
}

C++: #include <iostream>
#include <vector>
#include <algorithm>
#include <cmath>
#include <complex>
#include <utility>

using namespace std;

typedef complex<double> Point;

bool eqv(double a, double b) {
    return abs(a-b) < 1e-10;
}

bool inCircle(const Point &p, const Point &center, double radius) {
    double dist2 = norm(p-center);
    double rad2 = radius*radius;
    return dist2 < radius*radius || eqv(dist2, rad2);
}

int check(const Point &center, const vector<Point> &points) {
    const int N = points.size();
    int cnt = 0;
    for(int i = 0; i < N; ++i) {
        if(inCircle(points[i], center, 1)) ++cnt;
    }
    return cnt;
}

int test(int a, int b, const vector<Point> &points) {
    const Point v(points[a]-points[b]);
    const double hlen = abs(v)/2;

    if(hlen >= 1) return 0;
    const double plen = sqrt(1 - hlen*hlen);
    Point perp = v*Point(0, 1);
    perp /= abs(perp);
    const Point mid = (points[a]+points[b]) / 2.0;
    return max(check(mid + perp*plen, points), check(mid - perp*plen, points));
}

int main() {
    cin.tie(0);
    ios::sync_with_stdio(0);

    while(true) {
        int N;
        cin >> N;
        if(!N) break;

        vector<Point> points(N);
        for(int i = 0; i < N; ++i) {
            cin >> points[i].real() >> points[i].imag();
        }
        int ans = 0;
        for(int i = 0; i < N; ++i) {
            for(int j = 0; j < i; ++j) {
                ans = max(ans, test(i, j, points));
            }
        }
        if(N == 1) ans = 1;
        cout << ans << endl;
    }
    return 0;
}

C++: #include<bits/stdc++.h>
using namespace std;
const int N = 10;
typedef complex <double> point;
point a[N];
double kc[N][N];
int n;
int quetgoc(int i,double r,int n){
    vector<pair<double,bool> > goc;
    for(int j=0;j<n;j++){
        if(kc[i][j]<=2*r&&i!=j){
            double B=acos(kc[i][j]/(2*r));
            double A= arg(a[j]-a[i]);
            double alpha=A-B;
            double beta=A+B;
            goc.push_back(make_pair(alpha,true));
            goc.push_back(make_pair(beta,false));
        }
    }
    sort(goc.begin(),goc.end());
    int cnt=1,ans=1;
     vector<pair<double, bool> >::iterator it;
    for(it=goc.begin();it!=goc.end();++it){
        if((*it).second) cnt++;
        else cnt--;
        if(cnt>ans) ans=cnt;
    }
    return ans;
}
int maxpoint(point a[],int n,int r){
    for(int i=0;i<n-1;i++)
        for(int j=i+1;j<n;j++)
        kc[i][j]=kc[j][i]=abs(a[i]-a[j]);
    int ans=0;
    for(int i=0;i<n;i++)
        ans=max(ans,quetgoc(i,r,n));
    return ans;
}
signed main()
{
    double x[10000],y[10000];
    int t;
    while(cin>>t&&t!=0){
    for(int i=0;i<t;i++) a[i]=0;
    for(int i=0;i<t;i++) cin>>x[i]>>y[i];
    for(int i=0;i<t;i++) a[i]=point(x[i],y[i]);
    int r = 1;
    cout << maxpoint(a,t, r)<<endl;
    }
}

C++: #include <cstdio>
#include <vector>
#include <cmath>

using namespace std;

const double r = 1.0;
const double PI = acos(-1.0);

struct Point {
    double x, y;
    Point() {}
    Point(double x, double y) : x(x), y(y) {}
};
Point operator+(const Point &a, const Point &b) {
    return Point(a.x + b.x, a.y + b.y);
}
Point operator-(const Point &a, const Point &b) {
    return Point(a.x - b.x, a.y - b.y);
}
Point operator*(const Point &a, double b) {
    return Point(a.x * b, a.y * b);
}
inline double Distance2(const Point &a, const Point &b) {
    return (a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y);
}
// p1??¨??????p1p2???theta?????¢?????????????????£????????????????????????????????????????????????
Point Circle(const Point &p1, const Point &p2, double theta) {
    Point p12 = p2 - p1;
    Point rp(p12.x * cos(theta) - p12.y * sin(theta), p12.x * sin(theta) + p12.y * cos(theta));
    double len = sqrt(p12.x * p12.x + p12.y * p12.y);

    return (rp * (1.0 / len)) + p1;
}

int main()
{
    int n;

    while (true) {
        // Input
        scanf("%d", &n);
        if (n == 0) break;
        vector<Point> p(n);
        for (int i = 0; i < n; ++i)
            scanf("%lf %lf", &p[i].x, &p[i].y);

        // Solve
        // 2???????????????????§£???????£???¨??????
        int ans = 1;
        for (int i = 0; i < n; ++i)
            for (int j = i + 1; j < n; ++j) {
                double d = Distance2(p[i], p[j]);

                if (d <= 4.0 * r * r) {
                    double theta = (PI - acos((2.0 - d) * 0.5)) * 0.5;

                    for (int s = -1; s <= 1; s += 2) {
                        Point c = Circle(p[i], p[j], s * theta);

                        int num = 2;
                        for (int k = 0; k < n; ++k) {
                            if (k != i && k != j && Distance2(c, p[k]) <= 1.0)
                                ++num;
                        }

                        ans = max(ans, num);
                    }
                }
            }

        // Output
        printf("%d\n", ans);
    }

    return 0;
}