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 <cassert>
#include <algorithm>
#include <bitset>
#include <cctype>
#include <cmath>
#include <complex>
#include <chrono>
#include <deque>
#include <functional>
#include <iomanip>
#include <iostream>
#include <map>
#include <numeric>
#include <queue>
#include <set>
#include <stack>
#include <string>
#include <unordered_map>
#include <vector>
#include <random>
#include <memory>
#include <utility>
#include <limits>
#include "limits.h"
 
#define rep(i, a, b) for (long long (i) = (a); i < (b); i++)
#define all(i) i.begin(), i.end()
#define debug(...) std::cerr << "[" << #__VA_ARGS__ << "]:", debug_out(__VA_ARGS__)

void debug_out(){std::cerr<<std::endl;}

template<typename Head,typename... Tail>
void debug_out(Head h,Tail... t){
  std::cerr<<" "<<h;
  debug_out(t...);
}
 
template <typename T1, typename T2>
std::ostream& operator<<(std::ostream& os, std::pair<T1, T2> pa) {
  return os << pa.first << " " << pa.second;
}
 
template <typename T>
std::ostream& operator<<(std::ostream& os, std::vector<T> vec) {
  for (int i = 0; i < vec.size(); i++)os << vec[i] << (i + 1 == vec.size() ? "" : " ");
  return os;
}
 
template<typename T1,typename T2>
inline bool chmax(T1& a,T2 b){return a<b && (a=b,true);}
 
template<typename T1,typename T2>
inline bool chmin(T1& a,T2 b){return a>b && (a=b,true);}
 
long long pow_mod(long long a, long long b, long long mod=-1) {
  if(b==0)return 1;
  if ((a == 0)||(mod!=-1&&(a+mod)%mod==0)) {
    return 0;
  }
  long long x = 1;
  while (b > 0) {
    if (b & 1) {
      x = (mod!=-1)?(x * a) % mod:x*a;
    }
    a = (mod!=-1)?(a * a) % mod:a*a;
    b >>= 1;
  }
  return x;
}

// const long long MOD = 998244353;
const long long MOD = 1e9 + 7;

using ll = long long;
using P=std::pair<long long,long long>;

using Point=std::complex<long double>;

const long double EPS=1e-10;
const long double PI=std::acos(-1.0l);

bool EQ(long double a,long double b){return std::abs(a-b)<EPS;}
bool EQP(Point a,Point b){return EQ(a.real(),b.real())&&EQ(a.imag(),b.imag());}

Point unitvec(Point a){return a/std::abs(a);}

int main() {
  std::cin.tie(nullptr);
  std::ios::sync_with_stdio(false);

  ll n;
  while(std::cin>>n,n){
    std::vector<Point> point;
    rep(i,0,n){
      long double a,b;
      std::cin>>a>>b;
      point.emplace_back(a,b);
    }

    std::vector<Point> candi;
    rep(i,0,n){
      rep(j,i+1,n){
        if(std::abs(point[j]-point[i])<=2.0l+EPS){
          Point mid=(point[i]+point[j])/2.0l;

          Point unit=unitvec(point[i]-mid);
          unit*=std::sqrt(1.0l-std::norm(point[i]-mid));
          candi.push_back(unit*Point(0,1)+mid);
          candi.push_back(unit*Point(0,-1)+mid);
        }
      }
    }

    ll ans=1;
    for(auto e:candi){
      ll cnt=0;
      rep(i,0,n){
        if(std::abs(e-point[i])<=1.0l+EPS)cnt++;
      }
      chmax(ans,cnt);
    }
    std::cout<<ans<<"\n";
  }

  return 0;
}

C++: #include<bits/stdc++.h>
using namespace std;

double distance(pair<double,double> p1,pair<double,double> p2){
	return sqrt((p1.first-p2.first)*(p1.first-p2.first)+(p1.second-p2.second)*(p1.second-p2.second));
}

int main(){
	int n;
	while(cin>>n,n){
		vector<double> x(n);
		vector<double> y(n);
		int res=1;
		for(int i=0;i<n;i++) cin>>x[i]>>y[i];
		for(int i=0;i<n;i++){
			for(int j=i+1;j<n;j++){
				pair<double,double> midp=make_pair((x[i]+x[j])/2,(y[i]+y[j])/2);
				double a=-(x[i]-x[j]),b=(y[i]-y[j]);
				if(a<0.000001) a=0;
				if(b<0.000001) b=0;
				//a(x-mid[i].first)+b(y-mid.second)=0
					for(int k=0;k<2;k++){
						pair<double,double> lef=midp;
						pair<double,double> rig;
						if(k==0) rig=make_pair(-b*1e8+midp.first,a*1e8+midp.second);
						if(k==1) rig=make_pair(b*1e8+midp.first,-a*1e8+midp.second);
						while(distance(lef,rig)>=0.00000001){
								pair<double,double> mid=make_pair((lef.first+rig.first)/2,(lef.second+rig.second)/2);
								if(distance(mid,make_pair(x[i],y[i]))>1){
									rig=mid;
								}
								else{
									lef=mid;
								}

						}
						int count=0;
						for(int r=0;r<n;r++){
							if(distance(make_pair(x[r],y[r]),lef)<1+0.00001){
								count++;
							}
						}
						res=max(res,count);
					}
			}

		}
		cout<<res<<endl;

	}

	return 0;
}

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;

double Dis(double x1, double y1, double x2, double y2);
int main(void){
	while(1){
		int n;
		int ans = 1;
		vector< pair<double, double> > p;

		cin >> n;
		if(n == 0) break;

		for(int i=0; i<n; i++){
			pair<double, double> t;
			cin >> t.first >> t.second;
			p.push_back(t);
		}

		for(int i=0; i<n; i++){
			for(int j=i+1; j<n; j++){
				//中点
				double mpx = (p[i].first + p[j].first) / 2.0;
				double mpy = (p[i].second + p[j].second) / 2.0;
				//2点間の距離
				double d1 = Dis(p[i].first, p[i].second, p[j].first, p[j].second);
				if(d1 > 2.0) continue;
				//p[i]からp[j]へ向かうベクトル
				double v1x = (p[j].first - p[i].first) / d1;
				double v1y = (p[j].second - p[i].second) / d1;
				//中点から円の中心へ向かう単位ベクトル
				double v2x = -v1y;
				double v2y = v1x;
				//中心から円の中心へ向かうベクトルの長さ
				double d2 = sqrt(1.0 - (d1*d1/4.0));

				//円の中心の座標
				double cx = mpx + v2x * d2;
				double cy = mpy + v2y * d2;

				int cnt = 0;
				for(int k=0; k<n; k++){
					if(k == i || k == j || (cx-p[k].first)*(cx-p[k].first) + (cy-p[k].second)*(cy-p[k].second) <= 1.0)	cnt++;
				}
				if(cnt > ans){
					ans = cnt;
				}

				cx = mpx - v2x * d2;
				cy = mpy - v2y * d2;

				cnt = 0;
				for(int k=0; k<n; k++){
					if(k == i || k == j || (cx-p[k].first)*(cx-p[k].first) + (cy-p[k].second)*(cy-p[k].second) <= 1.0)	cnt++;
				}
				if(cnt > ans){
					ans = cnt;
				}
			}
		}

		cout << ans << endl;


	}
	
	return 0;
}

double Dis(double x1, double y1, double x2, double y2){
	return sqrt((x1 - x2)*(x1 - x2) + (y1 - y2)*(y1 - y2));
}