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 <vector>
#include <string>
#include <complex>
#include <algorithm>
using namespace std;

typedef complex<double> p;

int n;
p ps[300];

int main() {
	while (cin >> n , n != 0) {
		for (int i = 0; i < n; i++) {
			double x, y;
			cin >> x >> y;
			ps[i] = p(x, y);
		}

		if (n == 1) {
			cout << 1 << endl;
			continue;
		}

		int ans = 0;
		for (int i = 0; i < n; i++) {
			for (int j = 0; j < n; j++) {
				p v = ps[i] - ps[j], ve = v / abs(v);
				double l = sqrt(1 - norm(v) * 0.25);
				p vn1 = ve * p(0, 1) * l, vn2 = ve * p(0, -1) * l;
				p c1 = ps[j] + v * 0.5 + vn1, c2 = ps[j] + v * 0.5 + vn2;
				int count1 = 2, count2 = 2;
				for (int k = 0; k < n; k++) {
					if (k == i || k == j) continue;
					if (norm(c1 - ps[k]) <= 1.0) count1++;
					if (norm(c2 - ps[k]) <= 1.0) count2++;
				}
				ans = max(ans, count1);
				ans = max(ans, count2);
			}
		}

		cout << ans << endl;
	}
}

C++: #include <cstdio>
#include <complex>
#include <vector>
using namespace std;
typedef complex<double> P;
const double PI=3.141592653589793238462643383279;

bool inC(P p, P c, double r){return abs(p-c)<r+1e-9;}

P getE(P p){return p/abs(p);}

double sq(double x){ return  x*x;}

int main(){
    int n,i,j,k,ans,cnt;
    double x, y;
    while(scanf("%d",&n)){
        if(n==0) break;
        ans=0;
        vector<P> ps(n);
        for(i=0; i<n; i++){
            scanf("%lf%lf",&x,&y);
            ps[i] = P(x,y);
        }
        for(i=0;i<n;i++)
            for(j=i+1;j<n;j++)if(abs(ps[i]-ps[j])<2){
                double d=sqrt(1-sq(abs(ps[j]-ps[i])*0.5));
                P c1=0.5*(ps[i]+ps[j])+getE(getE(ps[j]-ps[i]))*polar(d,PI/2);
                P c2=0.5*(ps[i]+ps[j])+getE(getE(ps[j]-ps[i]))*polar(d,-PI/2);
                cnt=0;
                for(k=0;k<n;k++)if(inC(ps[k],c1,1.0))cnt++;
                if(cnt>ans)ans=cnt;
                cnt=0;
                for(k=0;k<n;k++)if(inC(ps[k],c2,1.0))cnt++;
                if(cnt>ans)ans=cnt;
            }
        printf("%d\n",ans);
    }
    return 0;
}

C++: #include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>
#include <queue>
#include <stack>
#include <set>

using namespace std;


int main(){
	int n;
	double x1,y1,xc,yc,dis,vx,vy,ex,ey;
	while(1){
		cin >> n;
		if(n==0)break;
		vector<double> x(n);
		vector<double> y(n);
		for(int i=0;i<n;i++){
			cin >> x[i] >> y[i];
		}
		if(n==1){
			cout << 1 << endl;
			continue;
		}
		int ans =0;
		for(int i=0;i<n;i++){
			for(int j=0;j<n;j++){
				if(i==j)break;
				xc = (x[i]+x[j])/2;
				yc = (y[i]+y[j])/2;
				if((xc-x[i])*(xc-x[i])+(yc-y[i])*(yc-y[i])>0.9999)continue;
				dis = sqrt(1.0-(xc-x[i])*(xc-x[i])-(yc-y[i])*(yc-y[i]));
				vx = x[i]-x[j];
				vy = y[i]-y[j];
				ex = dis * vy / sqrt(vx*vx+vy*vy);
				ey = -dis* vx / sqrt(vx*vx+vy*vy);
				x1 = xc + ex;
				y1 = yc + ey;
				int count=0;
				for(int k=0;k<n;k++){
					if(k==i||k==j)count +=1;
					else{
						if((x1-x[k])*(x1-x[k])+(y1-y[k])*(y1-y[k])<=1.0001){
							count +=1;
						}
					}
				}
				ans = max(ans,count);
				x1 = xc - ex;
				y1 = yc - ey;
				count = 0;
				for(int k=0;k<n;k++){
					if(k==i||k==j)count +=1;
					else{
						if((x1-x[k])*(x1-x[k])+(y1-y[k])*(y1-y[k])<=1.0001){
							count +=1;
						}
					}
				}
				ans = max(ans,count);
			}
		}
		cout << ans << endl;
	}
	return 0;
}

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

using namespace std;

int count( vector<complex<double>> &vpos );

int main(int argc, char const *argv[]) {
    int n, result;
    double x, y;
    vector<complex<double>> vpos;
    vector<int> results;

    while(1){
        cin>>n;
        if( n == 0 ) break;
        for(int i = 0; i < n; i++){
            cin>>x>>y;
            vpos.push_back( complex<double> (x,y) );
        }
        results.push_back( count( vpos ) );
    }

    for (auto it = results.begin(); it != results.end(); it++) {
        cout<<*it<<endl;
    }

    return 0;
}

int count( vector<complex<double>> &vpos ){
    int count = 0;
    int tmp_count;
    complex<double> m, n, center;

    for(auto fi = vpos.begin(); fi != vpos.end(); fi++){
        for(auto sc = vpos.begin(); sc != vpos.end(); sc++){
            if( (*fi).real() == (*sc).real() ) continue;
            tmp_count = 0;

            m = ( (*fi) + (*sc) ) / 2.0;
            n = ( (*fi) -  (*sc) ) * complex<double> ( 0.0, 1.0 );
            n = n / abs(n);

            double x = abs( (*fi) -  (*sc)  ) / 2.0;
            center =  m + n * sqrt( 1.0 - x * x );

            for(auto it = vpos.begin(); it != vpos.end(); it++){
                double dist = abs( center - (*it ) );
                if( dist <  2.0 ){
                    tmp_count++;
                }
            }
            if( count < tmp_count ){
                count = tmp_count;
            }
        }
    }
    return count;
}

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 long 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 = 0;
    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-EPS, 1));
                v.push_back(Pi(2*PI+EPS, -1));
                v.push_back(Pi(0-EPS, 1));
                v.push_back(Pi(r+EPS, -1));
            } else {
                v.push_back(Pi(l-EPS, 1));
                v.push_back(Pi(r+EPS, -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;
}