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 <vector>
#include <list>
#include <map>
#include <set>
#include <deque>
#include <stack>
#include <bitset>
#include <algorithm>
#include <functional>
#include <numeric>
#include <utility>
#include <sstream>
#include <queue>
#include <iostream>
#include <iomanip>
#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <cctype>
#include <string>
#include <cstring>
#include <ctime>
#include <fstream>
#include <stdio.h>
#include <complex>

using namespace std;

//conversion
//------------------------------------------
inline int toInt(string s) {int v; istringstream sin(s);sin>>v;return v;}
template<class T> inline string toString(T x) {ostringstream sout;sout<<x;return sout.str();}

//typedef
//------------------------------------------
typedef vector<int> VI;
typedef vector<VI> VVI;
typedef vector<string> VS;
typedef pair<int, int> PII;
typedef pair<int, PII> TIII;
typedef long long LL;
typedef vector<LL> VLL;

//container util

//------------------------------------------
#define ALL(a)  (a).begin(),(a).end()
#define RALL(a) (a).rbegin(), (a).rend()
#define PB push_back
#define MP make_pair
#define SZ(a) int((a).size())
#define EACH(i,c) for(typeof((c).begin()) i=(c).begin(); i!=(c).end(); ++i)
#define EXIST(s,e) ((s).find(e)!=(s).end())
#define SORT(c) sort((c).begin(),(c).end())
#define MT(a,b,c) MP(a, MP(b, c))

//repetition
//------------------------------------------
#define FOR(i,a,b) for(int i=(a);i<(b);++i)
#define REP(i,n)  FOR(i,0,n)

double x[300], y[300];
int n;

int count(double x0, double y0, int p, int q){
	int ret = 2;
	REP(i,n){
		if(i==p||i==q)continue;
		double dx = x0-x[i];
		double dy = y0-y[i];
		if(pow(dx, 2.0)+pow(dy, 2.0)<1.0)ret++;
	}
	return ret;
}

int main(){
	while(1){
		int ret = 1;
		cin>>n;
		if(n==0)break;
		REP(i,n)scanf("%lf %lf", &x[i],&y[i]);
		REP(i,n){
			FOR(j, i+1,n){
				double dx = x[j]-x[i];
				double dy = y[j]-y[i];
				double d = sqrt(pow(dx, 2.0)+pow(dy, 2.0));
				if(d < 2.0){
					double mx = (x[j]+x[i])*0.5;
					double my = (y[j]+y[i])*0.5;
					{
						double x1 = mx-dy/d*sqrt(1.0-d*d/4.0);
						double y1 = my+dx/d*sqrt(1.0-d*d/4.0);
						ret = max(ret, count(x1, y1, i, j));
					}
					{
						double x1 = mx+dy/d*sqrt(1.0-d*d/4.0);
						double y1 = my-dx/d*sqrt(1.0-d*d/4.0);
						ret = max(ret, count(x1, y1, i, j));
					}
				}
			}
		}
		cout << ret << endl;
	}
	return 0;
}

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

using namespace std;

double EPS=1e-10;

double add(double a,double b){
  if(abs(a+b)<EPS * (abs(a)+abs(b)))return 0;
  return a+b;
}

struct point{
  double x,y;
  point(){}
  point(double x,double y) : x(x) , y(y) {}

  point operator + (point p) {
    return point(add(x,p.x),add(y,p.y));
  }
  point operator - (point p) {
    return point(add(x,-p.x),add(y,-p.y));
  }
  point operator * (double d) {
    return point(x*d,y*d);
  }
  point operator / (double d) {
    return point(x/d,y/d);
  }
  bool operator == (const point &p)const{
    return abs(x-p.x)<EPS && abs(y-p.y)<EPS;
  }
};

struct circle{
  point p;
  double r;
  circle(){}
  circle(point p,double r):p(p),r(r){}
};

double dist(point a, point b){
  return sqrt(pow(a.x-b.x,2)+pow(a.y-b.y,2));
}

vector<point> circle_cross_point(circle a, circle b){
  vector<point>ret;

  double x=a.p.x-b.p.x,y=a.p.y-b.p.y;
  double L=x*x+y*y;

  if(L>pow(a.r+b.r,2))return ret;

  L=sqrt(L);
  double C=atan2(b.p.y-a.p.y,b.p.x-a.p.x);
  double alpha=acos((L*L+a.r*a.r-b.r*b.r)/(2*L*a.r));

  point r1,r2;
  r1.x=a.p.x+cos(C+alpha),r1.y=a.p.y+sin(C+alpha);
  r2.x=a.p.x+cos(C-alpha),r2.y=a.p.y+sin(C-alpha);

  ret.push_back(r1),ret.push_back(r2);
  return ret;
}

int main(void){
  int n;
  point a;
  vector<point>P;
  vector<circle>C;

  while(cin >> n,n){
    P.clear();
    C.clear();
    for(int i=0;i<n;i++){
      cin >> a.x >> a.y;
      P.push_back(a);
    }

    for(int i=0;i<n;i++){
      for(int j=i+1;j<n;j++){
	if(dist(P[i],P[j])>2.0+EPS)continue;
	vector<point>res=circle_cross_point(circle(P[i],1.0),circle(P[j],1.0));
	for(int k=0;k<res.size();k++)
	  C.push_back(circle(res[k],1));
      }
    }

    int ans=1;
    for(int i=0;i<C.size();i++){
      int cnt=0;
      for(int j=0;j<P.size();j++)
	if(dist(C[i].p,P[j])<1.0+EPS)cnt++;

      ans=max(ans,cnt);
    }
    cout << ans << endl;
  }

  return 0;
}

C++: #include <cstdio>
#include <cmath>
#include <cstring>
#include <cassert>
#include <iostream>
#include <string>
#include <vector>
#include <queue>
#include <stack>
#include <algorithm>
#include <map>

int n;
double x[512];
double y[512];

void init() {
  for(int i = 0; i < (int)512; ++i) {
    x[i] = y[i] = 0.0;
  }
} 

inline double Heron(double a, double b, double c) {
  double s = (a + b + c) / 2;
  return sqrt(s*(s-a)*(s-b)*(s-c));
}

class solve_ {
public:
  int operator() (void) {
    int res = 0;
    for(int i = 0; i < (int)n; ++i) {
      for(int j = 0; j < (int)i; ++j) {
        res = std::max(res, CountPoints(i, j));
      }
    }
    return res;
  }
  int CountPoints(int i, int j) {
    double dx = x[i] - x[j];
    double dy = y[i] - y[j];
    double dist = sqrt(dx * dx + dy * dy);
    auto t1 = std::pair<double, double>((x[i] + x[j]) / 2, (y[i] + y[j]) / 2);
    auto t2 = std::pair<double, double>(dx / dist, dy / dist);
    double h = 2 * Heron(1, 1, dist) / dist;
    auto t3 = std::pair<double, double>(t2.first * h, t2.second * h);
    auto c1 = std::pair<double, double>(t1.first + t3.second, t1.second - t3.first);
    auto c2 = std::pair<double, double>(t1.first - t3.second, t1.second + t3.first);
    int count1 = 0;
    int count2 = 0;
    {
      double cx = c1.first, cy = c1.second;
      for(int k = 0; k < (int)n; ++k) {
        if( k == i or k == j ) continue;
        if( (cx - x[k]) * (cx - x[k]) + (cy - y[k]) * (cy - y[k]) <= 1 ) {
          count1 += 1;
        }
      }
    }
    {
      double cx = c2.first, cy = c2.second;
      for(int k = 0; k < (int)n; ++k) {
        if( k == i or k == j ) continue;
        if( (cx - x[k]) * (cx - x[k]) + (cy - y[k]) * (cy - y[k]) <= 1 ) {
          count2 += 1;
        }
      }
    }
    int res = std::max(count1, count2) + 2;
    //printf("(%d, %d) -> %d\n", i, j, res);
    return res;
  }
};
solve_ solve;

int main() {
  for(;;) {
    init();
    scanf("%d", &n);
    if( n == 0 ) break;
    for(int i = 0; i < (int)n; ++i) {
      scanf("%lf %lf", &x[i], &y[i]);
    }
    printf("%d\n", solve());
  }
  return 0;
}