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<string>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<stack>
#include<queue>
#include<cmath>
#include<numeric>
#include<cassert>
#include<limits>
#include<iomanip>
using namespace std;
typedef unsigned long long int ull;
typedef vector<char> vc;
typedef vector<vc> vvc;



typedef pair<double,double> Point;
#define F first
#define S second
#define MK make_pair

Point pls(Point a,Point b){
    return MK(a.F+b.F,a.S+b.S);
}

Point mns(Point a,Point b){
    return MK(a.F-b.F,a.S-b.S);
}

Point T(Point t){
    return MK(t.S,-1.0*t.F);
}

double abs(Point a){
    return sqrt(a.F*a.F+a.S*a.S);
}

double norm(Point a){
    return a.F*a.F+a.S*a.S;
}

Point mid(Point a,Point b){
    return MK((a.F+b.F)/2,(a.S+b.S)/2);
}

bool in(Point c,Point t){
    if(norm(mns(t,c))<=1){return true;}
    return false;
}


int main(){
    while(1){
        int n;
        cin>>n;
        if(n==0){break;}
        vector<Point> vp;
        int mx=0;
        for(int i=0;i<n;i++){
            double a,b;
            cin>>a>>b;
            vp.push_back(MK(a,b));
        }
        sort(vp.begin(),vp.end());
        for(int i=0;i<vp.size();i++){
            for(int t=i+1;t<vp.size() && vp[t].F-vp[i].F<=2;t++){
                if(abs(mns(vp[t],vp[i]))>2){continue;}
                Point m=mid(vp[i],vp[t]);
                Point l=mns(vp[t],vp[i]);
                Point lt=T(l);
                double h=sqrt(1.0-abs(l)*abs(l)/4.0);
                lt.F=h*lt.F/abs(l);
                lt.S=h*lt.S/abs(l);
                Point c=pls(m,lt);
                int count=2;
                for(int k=0;k<vp.size();k++){
                    if(k!=i && k!=t && in(c,vp[k])){count++;}
                }
                mx=max(mx,count);
                c=mns(m,lt);
                count=2;
                for(int k=0;k<vp.size();k++){
                    if(k!=i && k!=t && in(c,vp[k])){count++;}
                }
                mx=max(mx,count);
            }
        }
        cout<<mx<<endl;
    }
    return 0;
}

C++: #include<iostream>
#include<math.h>
#include<vector>
using namespace std;

// 点
class Point {
  public:
    float x, y;
    Point(double ix = 0.0, double iy = 0.0) {
      set(ix, iy);
    }
    void set(double ix = 0.0, double iy = 0.0) {
      x = ix;
      y = iy;
    }
    double length() {
      return sqrt(x * x + y * y);
    }
    double distanceTo(Point p) {
      double dx = p.x - x, dy = p.y - y;
      return sqrt(dx * dx + dy * dy);
    }
};

int main() {
  vector<int> results;
  while (true) {
    // 入力
    int n;
    cin >> n;
    if (n == 0) break;
    vector<Point> p(n);
    for (int i=0; i<n; i++) {
      double x, y;
      cin >> x >> y;
      p[i].set(x, y);
    }
    
    // とりあえず格子状に調べてみる（）
    const float X = 100;
    Point pc;
    int max = 0;
    for (int xi=0; xi<=X; xi++) {
      for (int yi=0; yi<=X; yi++) {
	float x = xi * 10 / X, y = yi * 10 / X; // 0 ~ 10に補正
	// 距離が1以下のものをカウント
	pc.set(x, y);
	int count = 0;
	for (int i=0; i<n; i++) {
	  if (pc.distanceTo(p[i]) <= 1) count++;
	}
	if (count > max) max = count;
      }
    }
    results.push_back(max);
  }
  
  int size = results.size();
  for (int i=0; i<size; i++) {
    cout << results[i] << endl;
  }
  return 0;
}

C++: #include <bits/stdc++.h>

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

using namespace std;
namespace geom{
#define at(i) ((*this)[i])
#define pb push_back
#define X real()
#define Y imag()
#define SELF (*this)
typedef long double R;
typedef complex<R> P;
const R EPS = 1e-8;
const R PI=3.14159265358979323846264338327950288;
enum { TURE = 1, FALSE = 0, BORDER= -1};
inline int sig(const R&x){return (abs(x)<EPS?0:x>0?1:-1);}
inline int less(const R &x,const R &y){return sig(x-y)?x<y:BORDER;}
inline R inp(const P &a,const P &b){return (conj(a)*b).X;}
inline R outp(const P &a,const P & b){return (conj(a)*b).Y;}
inline R norm(const P &p){return p.X*p.X+p.Y*p.Y;}
inline P unit(const P & p){return p/abs(p);}
inline P proj(const P &s,const P &t){return t*inp(s,t)/norm(t);}
struct L :public vector<P>{
  L(const P &p1,const P &p2){this->pb(p1);this->pb(p2);}
  L(){}
  P dir()const {return at(1)-at(0);}
  int online(const P &p) const {return !sig(outp(p-at(0),dir()));}
};
struct S: public L {
    S(const P &p1, const P &p2): L(p1,p2){}
    S(){}
    int online(const P &p) const{
        if(!sig(norm( p-at(0) )) || !sig(norm(p - at(1))) ) return BORDER;
        return !sig(outp(p-at(0), dir())) && inp(p-at(0), dir()) > EPS && inp(p-at(1),-dir()) > -EPS;
    }
};
struct C :public P{
  C(){}
  C(const P &p,const R r):P(p),r(r){}
  R r;
  int inside(const P &p)const { return less(norm(p-SELF),r*r);}
};
//inline P proj(const P &s,const L &t){return t[0] + proj(s-t[0], t[1]-t[0]);}

inline int intersect (const C &a,const C &b){
  return less((a.r-b.r)*(a.r-b.r),norm(a-b))+less(norm(a-b),(a.r+b.r)*(a.r+b.r))-1;
}

inline S crosspoint(const C &c1,const C &c2){
  if(!intersect(c1,c2)) return S();
  R d=abs(c1-c2);
  R x=(c1.r*c1.r-c2.r*c2.r+d*d)/(2*d);
  R h=sqrt(max<R>(0.,c1.r*c1.r-x*x));
  P u=unit(c2-c1);
  return S(c1+u*x+u*P(0,-1)*h,c1+u*x+u*P(0,1)*h);
}

// inline S crosspoint(const C &c,const L &l){
//   R d2=dist2(l,c);
//   if(c.r*c.r+EPS < d2) return S();
//   P m= proj(c,l);
//   P u = unit(l[1]-l[0]);
//   R d=sqrt(max<R>(.0,c.r*c.r-d2));
//   return S(m+u*d,m-u*d);
// }

S circlePPR(const P& a, const P& b, R r){
    return crosspoint( C(a,r), C(b,r) );
}
}
using namespace geom;

int main(){
    int N;
    while(cin>>N, N){
        int ans = 1;
        P po[305];
        rep(i,N){
            R x, y;
            cin >> x >> y;
            po[i] = P(x, y);
        }
        rep(i,N) rep(j,i){
            //if( i==j ) continue;
            S s = circlePPR( po[i], po[j], 1.0 );
            //if( s.empty() ) continue;
            //for(auto t: s) cout << t << " "; cout << endl;
            for(auto k: s){
                int count = 0;
                rep(h,N){
                    if( abs(k-po[h]) < 1.0 + EPS ){
                        count++;
                    }
                }
                ans = max(ans, count);
            }
        }
        cout << ans << endl;
    }
    return 0;
}