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 <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 = 0;
        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;
}

C++: #include <bits/stdc++.h>
using namespace std;
#define int long long   // <-----!!!!!!!!!!!!!!!!!!!
 
#define rep(i,n) for (int i=0;i<(n);i++)
#define rep2(i,a,b) for (int i=(a);i<(b);i++)
#define rrep(i,n) for (int i=(n)-1;i>=0;i--)
#define rrep2(i,a,b) for (int i=(a)-1;i>=b;i--)
#define all(a) (a).begin(),(a).end()
#define rall(a) (a).rbegin(),(a).rend()
#define printV(_v) for(auto _x:_v){cout<<_x<<" ";}cout<<endl
#define printVS(vs) for(auto x : vs){cout << x << endl;}
#define printVV(_vv) for(auto _v:_vv){for(auto _x:_v){cout<<_x<<" ";}cout<<endl;}
#define printP(p) cout << p.first << " " << p.second << endl
#define printVP(vp) for(auto p : vp) printP(p);
 
typedef long long ll;
typedef pair<int, int> Pii;
typedef tuple<int, int, int> TUPLE;
typedef vector<int> vi;
typedef vector<vi> vvi;
typedef vector<vvi> vvvi;
typedef vector<Pii> vp;
const int inf = 1e9;
const int mod = 1e9 + 7;
 
// #define double long double
 
typedef complex<double> P;
typedef vector<P> G;
#define here(g, i) g[i]
#define next(g, i) g[(i + 1) % g.size()]
#define prev(g, i) g[(i - 1 + g.size()) % g.size()]
const double EPS = 1e-6;
const double INF = 1e12;
const double PI = acos(-1);
 
P readP() {
    double x, y;
    cin >> x >> y;
    return P(x, y);
}
 
// rot p around q by theta (counter-clockwise)
P rotP(P p, P q, double theta) {
    p -= q;
    double x = p.real(), y = p.imag();
    p = P(x * cos(theta) - y * sin(theta), x * sin(theta) + y * cos(theta));
    p += q;
    return p;
}
 
P center(P p, P q) {
    // P v = q - p;
    // double theta = asin(abs(v) / 2);
    // v *= (1 / abs(v));
    // P u = rotP(v, p, theta);
    // return p + u;
 
    // P v = q - p;
    // P h = v * P(0., 1.) / abs(v);
    // return p + (v / 2.) + h;
 
    P v = q - p;
    P h = P(v.imag(), -v.real());
    h /= abs(h);
    h *= sqrt(1 - norm(v / 2.));
    P c = p + (v / 2.) + h;
    return c;
}
 
signed main() {
    std::ios::sync_with_stdio(false);
    std::cin.tie(0);
 
    int n;
    while (cin >> n, n) {
        vector<P> p(n);
        rep(i, n) p[i] = readP();
 
        int ans = 1;
        rep(i, n) {
            rep2(j, i + 1, n) {
                if (i == j) continue;
                if (norm(p[i] - p[j]) > 4. + EPS) continue;
                rep(s, 2) {
                    P c = (s ? center(p[i], p[j]) : center(p[j], p[i]));
                    int cnt = 2;
                    rep(k, n) {
                        if (k == i || k == j) continue;
                        if (norm(p[k] - c) < 1. + EPS) {
                            cnt++;
                        }
                    }
                    ans = max(ans, cnt);
                }
            }
        }
 
        cout << ans << endl;
    }
}