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 <iostream>
#include <cstdio>
#include <cassert>
#include <cstring>
#include <vector>
#include <valarray>
#include <array>
#include <queue>
#include <set>
#include <unordered_set>
#include <map>
#include <unordered_map>
#include <algorithm>
#include <cmath>
#include <complex>
#include <random>

using namespace std;
typedef long long ll;
typedef unsigned long long ull;

typedef double R;
typedef complex<R> P;

const R EPS = 1e-10;
const R PI = acos((R)(-1));

R radNorP(R x) {
    return fmod(fmod(x, 2*PI)+2*PI, 2*PI);
}

const int MN = 330;

P p[MN];
typedef pair<R, int> Pi;

bool solve() {
    int n;
    cin >> n;
    if (!n) return false;

    for (int i = 0; i < n; i++) {
        R x, y;
        cin >> x >> y;
        p[i] = P(x, y);
    }

    int res = 1;
    for (int i = 0; i < n; i++) {
        vector<Pi> v;
        for (int j = 0; j < n; j++) {
            if (i == j) continue;
            if (2.0 < abs(p[j]-p[i])) continue;
            R th = acos(abs(p[j]-p[i])/2.0);
            R l = radNorP(arg(p[j]-p[i])-th);
            R r = radNorP(arg(p[j]-p[i])+th);
            if (r < l) {
                v.push_back(Pi(l, 1));
                v.push_back(Pi(2*PI, -1));
                v.push_back(Pi(0, 1));
                v.push_back(Pi(r, -1));
            } else {
                v.push_back(Pi(l, 1));
                v.push_back(Pi(r, -1));
            }

        }

        sort(v.begin(), v.end());
        int sm = 1;
        for (Pi p: v) {            
            sm += p.second;
            res = max(res, sm);
        }
    }
    cout << res << endl;
    return true;
}

int main() {
    while (solve()) {}
    return 0;
}

C++: #include <stdio.h>
#include <math.h>
#include "stdlib.h"

#define MAXN 300
#define eps 1e-8

typedef struct {
    double x;
    double y;
} vec2;

vec2 input[MAXN];
uint num_of_input;

uint init() {
    uint N;
    scanf("%u", &N);
    num_of_input = N;


    if (!N)
        return 0;

    for (uint i = 0; i < N; i++) {
        double x, y;
        scanf("%lf%lf", &x, &y);

        input[i].x = x;
        input[i].y = y;
    }


    return 1;


}

double dist(vec2 a, vec2 b) {
    return sqrt(
            (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)
    );
}

vec2 get_center(vec2 a, vec2 b){
    vec2 ret,mid={.x=(a.x+b.x)/2,.y=(a.y+b.y)/2};
    double angle=atan2(a.y-b.y,a.x-b.x);
    double distance_from_mid_to_center=sqrt(1.0-dist(a,b)*dist(a,b)/4.0);

    ret.y=mid.y-distance_from_mid_to_center*cos(angle);
    ret.x=mid.x+distance_from_mid_to_center*sin(angle);

    return ret;
}

int main() {

    while (init()) {
        uint ans = 1;

        for (uint i = 0; i < num_of_input;i++) {
            for (uint j = i + 1; j < num_of_input; j++) {
                vec2 center;
                uint cur_count = 0;
                if (dist(input[i], input[j]) > 2.0) continue;

                center = get_center(input[i], input[j]);

                for (uint k = 0; k < num_of_input; k++) {
                    if (dist(center, input[k]) < 1.0 + eps)
                        cur_count++;
                }
                ans = ans > cur_count ? ans : cur_count;
            }
        }

        printf("%u\n", ans);
    }

    return 0;
}

C++: /*
??????????????????????????????3??\??????3???????????¨?????????????????????????????????????????¨????????????
*/
#include <iostream>
#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <cstring>
#include <string>
#include <vector>
#include <utility>
#include <complex>
#include <set>
#include <map>
#include <queue>
#include <stack>
#include <deque>
#include <tuple>
#include <bitset>
#include <algorithm>
using namespace std;
typedef long double ld;
typedef long long ll;
typedef vector<int> vint;
typedef pair<int, int> pii;
typedef pair<ll, ll> pll;
typedef pair<double, double> pdd;
typedef complex<double> comd;
#define rep(i,n)	for(int i=0;i<n;i++)
#define srep(i,a,n)	for(int i=a;i<n;i++)
#define REP(i,n)	for(int i=0;i<=n;i++)
#define SREP(i,a,n)	for(int i=a;i<=n;i++)
#define rrep(i,n)	for(int i=n-1;i>=0;i--)
#define RREP(i,n)	for(int i=n;i>=0;i--)
#define all(a)	(a).begin(),(a).end()
#define mp(a,b)	make_pair(a,b)
#define mt	make_tuple
#define fst	first
#define scn second
#define bucnt(x)	__buildin__popcount(x)
#define debug(x)	cout<<"debug: "<<x<<endl

const ll inf = (ll)1e9;
const ll mod = (ll)1e9 + 7;
const ld eps = 1e-9;
const int dx[] = { 0,1,0,-1 };
const int dy[] = { 1,0,-1,0 };

double Dist2(pdd a, pdd b) {
	return (a.fst - b.fst)*(a.fst - b.fst) + (a.scn - b.scn)*(a.scn - b.scn);
}

vector<pdd> getCenter(pdd a, pdd b) {
	vector<pdd> ret;
	pdd c = mp((a.fst + b.fst) / 2.0, (a.scn + b.scn) / 2.0);
	pdd v = mp(a.scn - b.scn, b.fst - a.fst);
	double D = sqrt(v.fst*v.fst + v.scn*v.scn);
	v.fst /= D;	v.scn /= D;
	D = sqrt(1 - Dist2(a, b)/4.0);
	v.fst *= D;	v.scn *= D;
	ret.push_back(mp(c.fst + v.fst, c.scn + v.scn));
	ret.push_back(mp(c.fst - v.fst, c.scn - v.scn));
	return ret;
}

int main() {
	while (true) {
		int n;	cin >> n;
		if (n == 0)	break;
		vector<pdd> p(n,mp(0,0));
		rep(i, n) cin >> p[i].fst >> p[i].scn;
		int ret = 0;
		rep(i, n) {
			rep(j, i) {
				if (Dist2(p[i], p[j]) <= 4.0) {
					vector<pdd> c = getCenter(p[i], p[j]);
					int tmp = 0;
					rep(k, n) {
						if (Dist2(c[0], p[k]) <= 1.0 + eps)	tmp++;
					}
					ret = max(ret, tmp);
					tmp = 0;
					rep(k, n) {
						if (Dist2(c[1], p[k]) <= 1.0 + eps)	tmp++;
					}
					ret = max(ret, tmp);
				}
			}
		}
		cout << ret << endl;
	}
	return 0;
}