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 <string.h>
#include <algorithm>
#include <iostream>
#include <math.h>
#include <assert.h>
#include <vector>
#include <complex>
#include <queue>

using namespace std;
typedef long long ll;
typedef complex<double> Point;
static const double EPS = 1e-9;
static const double PI = acos(-1.0);

#define REP(i, n) for (int i = 0; i < (int)(n); i++)
#define FOR(i, s, n) for (int i = (s); i < (int)(n); i++)
#define FOREQ(i, s, n) for (int i = (s); i <= (int)(n); i++)
#define FORIT(it, c) for (__typeof((c).begin())it = (c).begin(); it != (c).end(); it++)
#define DEC(i, s) for (int i = (s); i >= 0; i--)

#define SIZE(v) (int)((v).size())
#define MEMSET(v, h) memset((v), h, sizeof(v))
#define FIND(m, w) ((m).find(w) != (m).end())

int n;
double R;
Point points[3000];

struct Node {
  Point p;
  double w;
  int bound;
  Node(Point p, double w, int bound) : p(p), w(w), bound(bound) {;}
  bool operator<(const Node &rhs) const {
    return bound < rhs.bound;
  }
};

int main() {
  R = 1.0;
  while (scanf("%d", &n), n) {
    priority_queue<Node> que;
    REP(i, n) {
      double x, y;
      scanf("%lf %lf", &x, &y);
      points[i] = Point(x, y);
      que.push(Node(points[i], R, 10000));
    }
    int ans = 1;
    while (!que.empty()) {
      Node node = que.top();
      que.pop();
      if (node.bound <= ans) { break; }
      int lans = 0;
      int nbound = 0;
      REP(i, n) {
        double dist = abs(points[i] - node.p);
        if (dist < R) {
          lans++;
          nbound++;
        } else if (dist < R + node.w * 1.41421356) {
          nbound++;
        }
      }
      ans = max(ans, lans);
      if (nbound <= ans) { continue; }
      const double dx[4] = { -1, -1, 1, 1 };
      const double dy[4] = { -1, 1, -1, 1 };
      double nw = node.w / 2.0;
      REP(i, 4) {
        double nx = node.p.real() + nw * dx[i];
        double ny = node.p.imag() + nw * dy[i];
        que.push(Node(Point(nx, ny), nw, nbound));
      }
    }
    printf("%d\n", ans);
  }
}

C++: #include <stdio.h>
#include <algorithm>
#include <iostream>
#include <stack>
#include <vector>
#include <math.h>
#include <queue>

using namespace std;

struct Point{
	double x,y;
};

void func(int N){

	Point point[N];
	double x_min = 100,y_min = 100,x_max = -1,y_max = -1;

	for(int i = 0; i < N; i++){
		scanf("%lf %lf",&point[i].x,&point[i].y);
		x_min = min(point[i].x,x_min);
		y_min = min(point[i].y,y_min);
		x_max = max(point[i].x,x_max);
		y_max = max(point[i].y,y_max);
	}

	if(x_min > 1)x_min -= 1;
	if(y_min > 1)y_min -= 1;
	if(x_max < 9)x_max += 1;
	if(y_max < 9)y_max += 1;

	int count,max_count = 0;

	for(double y = y_min; y <= y_max; y += 0.01){
		for(double x = x_min; x <= x_max; x += 0.01){
			count = 0;
			for(int i = 0; i < N; i++){
				if(sqrt((x-point[i].x)*(x-point[i].x)+(y-point[i].y)*(y-point[i].y)) <= 1.0)count++;
			}
			max_count = max(max_count,count);
		}
	}

	printf("%d\n",max_count);
}


int main(){

	int N;

	while(true){
		scanf("%d",&N);
		if(N == 0)break;

		func(N);
	}


    return 0;
}

C++: #include <cstdio>
#include <cstdlib>
#include <cmath>
#include <climits>
#include <cfloat>
#include <map>
#include <utility>
#include <set>
#include <iostream>
#include <memory>
#include <string>
#include <vector>
#include <algorithm>
#include <functional>
#include <sstream>
#include <complex>
#include <stack>
#include <queue>
#include <cstring>
#include <sstream>
#include <cassert>
#include <list>
using namespace std;
static const double EPS = 1e-5;
typedef long long ll;
typedef pair<int,int> PI;
typedef vector<int> vi;
#ifndef M_PI
#define M_PI 3.14159265358979323846
#endif
#define rep(i,n) for(int i=0;i<(int)n;++i)
#define FOR(i,c) for(__typeof((c).begin())i=(c).begin();i!=(c).end();++i)
#define ALL(c) (c).begin(), (c).end()
#define mp(a,b) make_pair(a,b)
#define pb(a) push_back(a)
#define F first
#define S second


int dx[]={0,1,0,-1};
int dy[]={1,0,-1,0};


main(){
  int n;
  while(cin>>n,n){
    double x[n],y[n];

    rep(i,n)cin>>x[i]>>y[i];

    int ans=0;

    rep(i,n){
      rep(j,i){
	double sx,sy;
	double tx=x[i]-x[j],ty=y[i]-y[j];
	double len=sqrt(tx*tx+ty*ty);
	if(len>2)continue;
	len=sqrt(1-len*len/4);
	double theta=atan2(ty,tx)+M_PI/2;
	int tans=0;
	sx=(x[i]+x[j])/2+cos(theta)*len;
	sy=(y[i]+y[j])/2+sin(theta)*len;
	rep(k,n){
	  tx=sx-x[k];
	  ty=sy-y[k];
	  if(tx*tx+ty*ty<=1+EPS)++tans;
	}
	ans=max(ans,tans);
	tans=0;
	sx=(x[i]+x[j])/2-cos(theta)*len;
	sy=(y[i]+y[j])/2-sin(theta)*len;
	rep(k,n){
	  tx=sx-x[k];
	  ty=sy-y[k];
	  if(tx*tx+ty*ty<=1+EPS)++tans;
	}
	ans=max(ans,tans);
      }
    }
    cout<<ans<<endl;
  }
}

C++: #include <iostream>
#include <math.h>

using namespace std;

int n;

double xi[300];
double yi[300];

double distant2(double x1,double y1,double x2,double y2)
{
	return (x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2);
}

int search(int i,int j)
{
	double mx = (xi[i] + xi[j]) / 2,my = (yi[i] + yi[j]) / 2;
	double diff2 = distant2(mx,my,xi[i],yi[i]);

	if(diff2 > 1.0)
		return 1; //1ツ古つつセツつッ

	double len = sqrt(1 - diff2);
	double diff = sqrt(diff2);
	double dx = len * (xi[i] - mx) / diff,dy = len * (yi[i] - my) / diff;

	int sign[] = {1,-1};
	int counter[2] = {0};
	for(int s = 0; s < 2; s++)
	{
		double x = mx + sign[s] * dx,y = my + sign[s] * dy;
		for(int k = 0; k < n; k++){
			if(k == i || k == j){counter[s]++; continue;}
			if(distant2(x,y,xi[k],yi[k]) < 1.0)
				counter[s]++;
		}
	}

	return max(counter[0],counter[1]);
}



int solve()
{
	for(int i = 0; i < n; i++){
		cin >> xi[i] >> yi[i];
	}

	int count = 0;

	for(int i = 0; i < n; i++){
		for(int j = i+1; j < n; j++){
			count = max(count,search(i,j));
		}
	}

	return count;
}

int main() {
  while (cin >>n,n) {
    cout << solve() << endl;
  }

  return 0;
}