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];
				p ve = v / abs(v);
				double l = sqrt(1 - norm(v) / 4);
				p vn1 = ve * p(0,1) * l, vn2 = ve * p(0,-1) * l;
				p c1 = ps[j] + v / 2.0 + vn1, c2 = ps[j] + v / 2.0 + 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) count1++;
					if (norm(c2 - ps[k]) <= 1) count2++;
				}
				ans = max(ans, count1);
				ans = max(ans, count2);
			}
		}

		cout << ans << endl;
	}
}

C++: #include<bits/stdc++.h>
using namespace std;

#define int long long
typedef vector<int>vint;
typedef pair<int,int>pint;
typedef vector<pint>vpint;
#define rep(i,n) for(int i=0;i<(n);i++)
#define reps(i,f,n) for(int i=(f);i<(n);i++)
#define all(v) (v).begin(),(v).end()
#define each(it,v) for(__typeof((v).begin()) it=(v).begin();it!=(v).end();it++)
#define pb push_back
#define fi first
#define se second
template<typename A,typename B>inline void chmin(A &a,B b){if(a>b)a=b;}
template<typename A,typename B>inline void chmax(A &a,B b){if(a<b)a=b;}

const double PI=acos(-1);
const double EPS=1e-6;
typedef complex<double>Point;
Point in(){
    double x,y;
    scanf("%lf%lf",&x,&y);
    return Point(x,y);
}

int N;
Point P[333];

signed main(){

    Point r=polar(1.0,PI/2);
    while(scanf("%lld",&N),N){
        rep(i,N)P[i]=in();

        int ma=1;
        rep(i,N)rep(j,N){
            Point b=(P[i]+P[j])*0.5;
            Point v=(P[i]-P[j])*0.5*r;
            for(double d:{-1.0,1.0}){
                Point u=v*d*sqrt(1.0-abs(v)*abs(v));
                Point p=b+u;
                int cnt=0;
                rep(k,N)if(abs(P[k]-p)<1.0+EPS)cnt++;
                chmax(ma,cnt);
            }
        }
        printf("%lld\n",ma);
    }
    return 0;
}

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

using namespace std;

typedef complex<double> xy_t;

const xy_t INF_VECT = xy_t(100, 100);
const xy_t ZERO_VECT = xy_t(0, 0);
const double ERR = 0.00001;

xy_t center(xy_t p, xy_t q);

int main()
{
//    ifstream fin("1132_input.txt");

    int num;
    while(1){
        cin >> num;
//        fin >> num;
        if(num == 0) break;
        
        xy_t point[num];
        for(int i = 0; i < num; i++){
            double x, y;
            cin >> x >> y;
//            fin >> x >> y;
            point[i] = xy_t(x, y);
        }

        if(num == 1){
            cout << "1" << endl;
        }

        else if(num > 1){
            int max = 1;
            for(int i = 0; i < num-1; i++){
                for(int j = i+1; j < num; j++){
                    xy_t vec = center(point[i], point[j]);
                    if(vec != INF_VECT){
                        xy_t mid = (point[i] + point[j]) * xy_t(0.5, 0);
                        xy_t cent = mid + vec;
                        
                        int cnt = 0;
                        for(int k = 0; k < num; k++){
                            xy_t r = point[k] - cent;
                            if(abs(r) <= 1 + ERR) cnt++;
                        }
                        if(cnt > max) max = cnt;

                        cent = mid - vec;
                        cnt = 0;
                        for(int k = 0; k < num; k++){
                            xy_t r = point[k] - cent;
                            if(abs(r) <= 1 + ERR) cnt++;
                        }
                        if(cnt > max) max = cnt;
                    }
                }
            }
            cout << max << endl;
        }
    }
    return 0;
}


xy_t center(xy_t p, xy_t q)
{    
    xy_t r = q-p;
    double d = abs(r);

    if(d > 2){
        return INF_VECT;
    }
    else if(d == 2){
        return ZERO_VECT;
    }
    else{
        double l = sqrt(4 - d*d) / 2;
        xy_t r_ = r * xy_t(0, 1);
        r_ *= xy_t(l/d, 0);
        return r_;
   }
}

C++: #include <iostream>
#include <vector>
#include <algorithm>

using namespace std;

#define rep(i,n) for(int i = 0;i < n;i++)

#include <complex>
#include <math.h>

typedef double D;
typedef complex<D> P;
typedef pair<P, P> L;
typedef vector<P> VP;

const D EPS = 1e-9;
const D INF = 1e9;

#define X real()
#define Y imag()

#define LE(n,m)	((n) < (m) + EPS)
#define GE(n,m)	((n) + EPS > (m))
#define EQ(n,m) (abs((n)-(m)) < EPS)

VP cpr(P a, P b, D r) {
	VP cs;
	P h = (b - a) * 0.5;
	D d = abs(h);
	if (d == 0 || d > r) return cs;
	D N = sqrt(r*r - d*d);
	P n = h * P(0, 1) * (N / d);
	cs.push_back(a + h + n);
	if (N > 0)cs.push_back(a + h - n);
	return cs;
}

int main() {
	int n;
	while (cin >> n && n) {
		VP v(n);
		rep(i, n) {
			D x, y;
			cin >> x >> y;
			v[i] = P(x, y);
		}
		int ans = 1;
		rep(i, n)rep(j, n)if (i != j) {
			auto t = cpr(v[i], v[j], 1.0);
			for (auto p : t) {
				int count = 0;
				rep(k, n) {
					if (abs(p - v[k]) <= 1.0 + EPS) {
						count++;
					}
				}
				ans = max(ans, count);
			}
		}
		cout << ans << endl;
	}

	return 0;
}

C++: #include <iostream>
#include <cmath>
#define EPS 1.0e-5
using namespace std;

typedef struct {
  double x,y;
} P;

double dis(P p ){
  return sqrt( p.x*p.x + p.y*p.y );
}
double dis(P p1,P p2 ){
  P p;
  p.x = p1.x-p2.x;
  p.y = p1.y-p2.y;
  return dis( p );
}


int main(){
  int n, mx;
  P p[301], m, p1,p2, mv,dmv, c[2];
  double x;

  while( cin>>n && n ){
    for( int i=0;i<n;i++ )
      cin >> p[i].x >> p[i].y ;

    mx=0;

    for( int i=0;i<n;i++ ){
      for( int j=i+1;j<n;j++ ){
	if( dis(p[i],p[j])<2.0 ){
	  p1=p[i]; p2=p[j];
	  m.x = (p1.x+p2.x)/2.0;
	  m.y = (p1.y+p2.y)/2.0;

	  x = sqrt( 1-dis(p1,m) );
	  // ツ陳心ツ点 cツづーツ仰づ淞づゥ
	  if( fabs( p1.x-p2.x ) < EPS ){
	    mv.x=1.0;  mv.y=0.0;
	  }else{
	    double a1 = ( p2.y-p1.y )/( p2.x-p1.x );
	    mv.x=1.0; mv.y=-1.0/a1;
	  }
	  double scale = dis( mv );
	  dmv.x = x*mv.x / scale;
	  dmv.y = x*mv.y / scale;

	  c[0].x=m.x+dmv.x; c[0].y=m.y+dmv.y;
	  c[1].x=m.x-dmv.x; c[2].y=m.y-dmv.y;

	  for( int cn=0;cn<2;cn++ ){
	    int sum=0;
	    for( int k=0;k<n;k++ ){
	      if( dis( c[cn],p[k] ) < 1.0 )
		sum++;
	    }
	    if( mx<sum ) mx=sum;
	  }
	}
      }
    }
    cout << mx << endl;
  }

  return 0;
}