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<cassert>
#include<string>
#include<algorithm>
#include<cstdio>
#include<cstdlib>
#include<cmath>
#include<stack>
#include<queue>
#include<set>
#include<map>
#include<tuple>
#include<numeric>
using namespace std;

using pii = pair<int,int>;
using ll = long long;
#define rep(i, j) for(int i=0; i < (int)(j); i++)
#define repeat(i, j, k) for(int i = (j); i < (int)(k); i++)
#define all(v) v.begin(),v.end()
#define debug(x) cerr << #x << " : " << x << endl

const int INF = 1 << 31;
const ll INF_ll = 1LL << 62;
const double EPS = 1e-7;

// vector
template<class T> ostream& operator << (ostream &os , const vector<T> &v) {
    for(const T &t : v) os << "\t" << t; return os << endl;
}
template<class T> istream& operator >> (istream &is , vector<T> &v) {
    for(T &a : v) is >> a; return is;
}
// pair
template<class T, class U> ostream& operator << (ostream &os , const pair<T, U> &v) {
    return os << "<" << v.first << " " << v.second << ">";
}

using Coordinate = double;

bool eq(Coordinate a, Coordinate b) {
    return abs(a - b) < EPS;
}

/*
 * Point
 */
struct Point {
    Coordinate x,y;
    Point():x(0),y(0){}
    Point(Coordinate x,Coordinate y):x(x),y(y){}
    pair<Coordinate,Coordinate> to_pair() const { return make_pair(x,y); }
    bool operator == (const Point &a) const { return eq(x, a.x) && eq(y, a.y); }
    bool operator != (const Point &a) const { return !(*this == a); }
    
    Point operator + (const Point& p) const { return Point(x + p.x, y + p.y); }
    Point operator - (const Point& p) const { return Point(x - p.x, y - p.y); }
    Point operator - () const { return Point(-x, -y); }
    Point& operator += (const Point& p) { x += p.x; y += p.y; return *this; }
    Point& operator -= (const Point& p) { x -= p.x; y -= p.y; return *this; }
    // ???????????¨ x, y???
    bool operator <  (const Point &a) const { return x < a.x && y < a.y; }
    bool operator <= (const Point &a) const { return *this == a || *this < a;}
    bool operator >  (const Point &a) const { return !(*this <= a); }
    bool operator >= (const Point &a) const { return !(*this < a); }
};
template<typename T> Point operator * (T k, const Point &p) { return Point(k * p.x, k * p.y); }
template<typename T> Point operator / (const Point &p, T k) { return Point(p.x / k, p.y / k); }
double dot(const Point &a, const Point &b) { return a.x * b.x + a.y * b.y; }
double cross(const Point &a, const Point &b) { return a.x * b.y - a.y * b.x; }
double norm(const Point& a) { return sqrt( dot(a, a) ); }

istream &operator >> (istream& is, Point& p){ return is >> p.x >> p.y;}
ostream &operator << (ostream& os, Point& p){ return os << "(" << p.x << "," << p.y << ")";}

double distance(const Point& a, const Point& b) { return norm(a - b); }
double distance2(const Point& a, const Point& b) { return dot(a - b, a - b); }

Point rotate(Point &p, double degree) {
    double rad = degree * M_PI / 180;
    return Point(cos(rad) * p.x - sin(rad) * p.y,
                 sin(rad) * p.x + cos(rad) * p.y);
}
// ??£?????´???????????????
Point orthonormal(Point &p) {
    Point p2 = Point(-1 * p.y, p.x); //rotate(p, 90);
    return p2 / norm(p2);
}


enum class ccw_t {
    COUNTER_CLOCK_WISE = 1,
    CLOCK_WISE = -1,
    STRAIGHT_C_A_B = 2,
    STRAIGHT_A_B_C = -2,
    STRAIGHT_A_C_B = 0
};


ccw_t ccw(Point a,Point b,Point c){
    Point ab = b - a;
    Point ac = c - a;
    if( cross(ab,ac) > 0 ) return ccw_t::COUNTER_CLOCK_WISE;//a-b-c ???????¨???????
    if( cross(ab,ac) < 0 ) return ccw_t::CLOCK_WISE;        //a-b-c ????¨???????
    if( dot(ab,ac) < 0 ) return ccw_t::STRAIGHT_C_A_B;      //c-a-b
    if( norm(ab) < norm(ac) )return ccw_t::STRAIGHT_A_B_C;  //a-b-c or a==b
    return ccw_t::STRAIGHT_A_C_B;                           //a-c-b or b==c or a==c
}

/*
 * Segment
 */

struct Segment{
    Point a, b;
    Segment(){}
    Segment(const Point& a, const Point& b):a(a), b(b){}
    pair<Point,Point> to_pair() const { return make_pair(a,b); }
    bool operator == (const Segment& s) const { return to_pair() == s.to_pair();}
    bool operator != (const Segment& s) const { return to_pair() != s.to_pair();}
};

double length(const Segment &s){ return distance(s.a, s.b); }
istream &operator >> (istream& is, Segment& s){ return is >> s.a >> s.b;}
ostream &operator << (ostream& os, Segment& s){ return os << s.a << "->" << s.b;}

bool does_intersect(const Segment &a, const Segment &b){
    return (int)ccw(a.a, a.b, b.a) * (int)ccw(a.a, a.b, b.b) <= 0 &&
            (int)ccw(b.a, b.b, a.a) * (int)ccw(b.a, b.b, a.b) <= 0;
}

double distance( const Segment &s , const Point &p){
    if( dot( s.b - s.a , p - s.a ) < EPS ) return norm( p - s.a );
    if( dot( s.a - s.b , p - s.b ) < EPS ) return norm( p - s.b );
    return abs(cross( s.b - s.a , p - s.a ) / norm( s.b - s.a ));
}

Point midPoint(const Segment &s) {
    return Point((s.a.x + s.b.x) / 2, (s.a.y + s.b.y) / 2);
}


bool solve() {
    int n; cin >> n;
    if(n == 0) return 0;
    vector<Point> P(n); cin >> P;
    int ans = 0;
    rep(i, n) rep(j, i) if(distance(P[i], P[j]) <= 2.0 + EPS){
        Point a = P[i] - P[j];
        // double len = sqrt(1 - pow(distance(P[i], P[j]) / 2, 2));
        double len = sqrt(1 - distance2(midPoint(Segment(P[i], P[j])), P[j]));
        Point c = len * orthonormal(a) + midPoint(Segment(P[i], P[j]));
        int cnt = 0;
        rep(l, n) {
            cnt += distance2(P[l], c) <= 1.0 + EPS;
        }
        ans = max(ans, cnt);
    }
    cout << ans << endl;
    return 1;
}

int main() {
    cin.tie(0);
    ios::sync_with_stdio(false);
    while(solve());
    return 0;
}