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++: // http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=1132&lang=jp

#include <iostream>
#include <complex>
#include <cmath>
#include <iomanip>
using namespace std;

typedef complex<double> xy_t;

const xy_t rot[2] = { xy_t(0, 1), xy_t(0, -1) };
const double eps = 1e-6; // 最初から誤差の見積もりがあるのでそれに甘えた。

int main() {
  int N;
  while (scanf("%d", &N) != EOF && N) {
    xy_t P[310];
    for (int i=0; i<N; i++) {
      double X, Y;
      scanf("%lf%lf", &X, &Y);
      P[i] = xy_t(X, Y);
    }
    int ans = 1;
    for (int i=0; i<N; i++) {
      for (int j=i+1; j<N; j++) {
	if ( abs(P[i] - P[j]) > 2 - eps ) continue;
	for (int k=0; k<2; k++) {
	  xy_t c = (P[i] + P[j])/2.0; // 中点
	  xy_t n = (P[i] - P[j]);
	  xy_t m = n/abs(n) * rot[k] * sqrt(1 - abs(n/2.0) * abs(n/2.0));
	  xy_t center = c + m;
	  int count = 0;
	  for (int k=0; k<N; k++) {
	    if (abs(center - P[k]) < 1 + eps) count++;
	    // cerr << fixed << setprecision(10) << abs(center - P[k]) << endl;
	  }
	  ans = max(ans, count);
	}
      }
    }
    printf("%d\n", ans);
  }
}

C++: #include <iostream>
#include <complex>
#include <cmath>
#include <algorithm>
using namespace std;
typedef complex<double> C;
#define EPS (1e-6)

int n;
C c[301];

void solve() {
  for(int i=0;i<n;i++)
  cin>>c[i].real()>>c[i].imag();
  
  int res = 1;
  for (int i = 0; i < n; i += 1) {
    for (int j = i+1; j < n; j += 1) {
      C v = c[j] - c[i];
      double absv = abs(v);
      if (absv > 2.0 + EPS) continue;
      // _
      C m = c[i] + (v / 2.0);
      // ~ÜÅÌ£
      double dist = sqrt(1.0 - absv * absv / 4.0);
      // PÊ@üxNg
      C unv  = (v * C(0, 1)) / absv;
      // S
      C nv = unv * dist;
      for (int k = 0; k < 2; k += 1) {          
        int cnt = 0;
        nv *= -1;
        C tc = m + nv;
        for (int l = 0; l < n; l += 1) {
          if (abs(tc - c[l]) < 1.0 + EPS) cnt++;
        }
        res = max(res, cnt);
      }
    }
  }
  cout<<res<<endl;
}

int main (int argc, char const* argv[]) {
  while (cin>>n, n) {
    solve();
  }
  return 0;
}

C++: #include <bits/stdc++.h>
using namespace std;
/*{{{*/  //template
#define rep(i,n) for(int i=0;i<n;i++)
#define INF 1<<29
#define LINF LLONG_MAX/3
#define mp make_pair
#define pb push_back
#define EB emplace_back
#define fi first
#define se second
#define all(v) ALL(v)
#define sz(x) (int)(x).size()
#define debug(x) cerr<<#x<<":"<<x<<endl
#define debug2(x,y) cerr<<#x<<","<<#y":"<<x<<","<<y<<endl
//struct fin{ fin(){ cin.tie(0); ios::sync_with_stdio(false); } } fin_;
struct Double{ double d; explicit Double(double x) : d(x){} };
ostream& operator<<(ostream& os,const Double x){ os << fixed << setprecision(20) << x.d; return os; }
template<typename T> ostream& operator<<(ostream& os,const vector<T>& vec){ os << "["; for(const auto& v : vec){ os << v << ","; } os << "]"; return os; }
template<typename T,typename U> ostream& operator<<(ostream& os,const pair<T,U>& p){ os << "(" << p.first << ","<< p.second <<")"; return os; }
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
typedef vector<int> vi;
typedef vector<vi> vvi;
ll gcd(ll a,ll b){ if(b==0) return a; else return gcd(b,a%b); }
constexpr double eps = 1e-14; 
constexpr ll mod = 1e9+7;
const int dx[]={1,0,-1,0} ,dy[] = {0,1,0,-1};
/*}}}*/

int N;

void solve(){
    vector<double> x(N),y(N);
    rep(i,N){ cin >> x[i] >> y[i]; }

    int ans=1;
    rep(i,N) rep(j,N) if(i!=j){
        //cout << "i,j = " << i << "," << j << endl;
        double vx = x[i]-x[j];
        double vy = y[i]-y[j];
        double d = sqrt(vx*vx+vy*vy);
        if(d+eps>2) continue;
        double r = sqrt(1-(d/2)*(d/2));

        double rx =  -vy*r/d;
        double ry = vx*r/d;

        double cx = (x[i]+x[j])/2 + rx;
        double cy = (y[i]+y[j])/2 + ry;

        int cnt=0;
        rep(k,N){
            double tx = cx-x[k];
            double ty = cy-y[k];
            double dd = sqrt(tx*tx+ty*ty);
            //cout << "dd=" << dd << endl;
            if(dd<=1+eps) cnt++;
        }
        //cout << "cnt : " << cnt << endl;
        ans = max(ans,cnt);
    }
    cout << ans << endl;
}

int main(){
    while(cin>>N){
        if(N==0) break;
        solve();
    }
}

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

using namespace std;

#define FOR(i,k,n) for(int i=(k); i<(int)(n); ++i)
#define REP(i,n) FOR(i,0,n)
#define FORIT(i,c) for(__typeof((c).begin())i=(c).begin();i!=(c).end();++i)

template<class T> void debug(T begin, T end){ for(T i = begin; i != end; ++i) cerr<<*i<<" "; cerr<<endl; }
inline bool valid(int x, int y, int W, int H){ return (x >= 0 && y >= 0 && x < W && y < H); }

typedef long long ll;
const int INF = 100000000;
const double EPS = 1e-8;
const int MOD = 1000000007;
int dx[8] = {1, 0, -1, 0, 1, -1, -1, 1};
int dy[8] = {0, 1, 0, -1, 1, 1, -1, -1};
#include <complex>
typedef complex<double> Point;
int main(){
  int N;
  while(cin >> N && N){
    vector<Point> ps(N);
    REP(i, N) cin >> ps[i].real() >> ps[i].imag();
    int ans = 0;
    REP(i, N)REP(j, N) if(i != j){
      double d = abs(ps[i] - ps[j]) / 2;
      if(d > 1) continue;
      double t = sqrt(1 - d * d);
      Point start = ps[i] + (ps[j] - ps[i]) * 0.5;
      Point vect = (ps[j] - ps[i]) * Point(0.0, 1.0);
      vect /= abs(vect);
      for(int sign = -1; sign <= 1; sign += 2){
        Point center = start + vect * t * (double)sign;
        int sub_ans = 0;
        REP(k, N) if(abs(center - ps[k]) < 1.0 + EPS) sub_ans ++;
        ans = max(ans, sub_ans);
      }
    }
    cout << ans << endl;
  }
  return 0;
}