AOJ 2167 Find the Point
問題
n本の直線が与えられる。
全ての直線から等距離にある点を求めよ。
複数あるときはManyと、存在しないときはNoneと出力せよ。
誤差は10^-4の絶対誤差が許容される。
制約条件
n≦100
方針
変な場合わけしてやった。
まずは同一の直線をマージ。
直線が3本以下なら絶対にmany
平行でない3直線が存在したら、
- 一点で交わるとき、その点のみが候補
- 一点で交わらないとき、その3点の内心または傍心が候補
平行でない3直線が存在しなかったら、直線は絶対に4本以下。
二組の平行な直線からなっている。
なので候補は1通りに定まる。
候補を見つけたら、候補から全ての直線への距離が等しいかどうかを見る。
EPS=10^-8だと通らなくて、もう少しゆるくしたら通った。
二等分線を全て列挙して、二等分線同士の交点を見るO(n^4)の解法だとずいぶん実装が楽っぽかったのでそっちのほうがよかった。
ソースコード
const double INF = 1e12, EPS = 5e-7; typedef complex<double> P; namespace std{ bool operator<(const P& a, const P& b){ return real(a) != real(b) ? real(a) < real(b) : imag(a) < imag(b); } } inline double cross(const P& a, const P& b){ return imag(conj(a)*b); } inline double dot(const P& a, const P& b){ return real(conj(a)*b); } struct L : public vector<P>{ L(const P &a, const P &b) { push_back(a); push_back(b); } }; typedef vector<P> G; inline int ccw(P a, P b, P c) { b -= a; c -= a; if(cross(b, c) > 0) return +1; // counter clockwise if(cross(b, c) < 0) return -1; // clockwise if(dot(b, c) < 0) return +2; // c--a--b on line if(norm(b) < norm(c)) return -2; // a--b--c on line return 0; } P crosspoint(const L &l, const L &m) { double A = cross(l[1] - l[0], m[1] - m[0]); double B = cross(l[1] - l[0], l[1] - m[0]); if (abs(A) < EPS && abs(B) < EPS) return P(rand(), rand()); // same line if (abs(A) < EPS) return P(rand(), rand());// !!!PRECONDITION NOT SATISFIED!!! return m[0] + B / A * (m[1] - m[0]); } P projection(const L &l, const P &p) { double t = dot(p-l[0], l[0]-l[1]) / norm(l[0]-l[1]); return l[0] + t*(l[0]-l[1]); } double distanceLP(const L &l, const P &p) { return abs(p - projection(l, p)); } P incenter(const P &x, const P &y, const P &z){ P a = y - x, b = z - x; P d1 = a / abs(a) + b / abs(b); P d2 = -a / abs(a) + (b - a) / abs(b - a); return crosspoint(L(P(0, 0), d1), L(a, a + d2)) + x; } G excenter(const G &t){ G res; rep(i, 3){ P a = t[i], b = t[(i + 1) % 3], c = t[(i + 2) % 3]; res.pb(crosspoint(L(a, a + (a - c) / abs(a - c) + (b - a) / abs(b - a)), L(b, b + (b - c) / abs(b - c) + (a - b) / abs(a - b)))); } return res; } bool parallel(const L &l, const L &m){ return abs(cross(l[0] - l[1], m[0] - m[1])) < EPS; } int n; int main(){ while(cin >> n && n){ vector<L> ls; rep(i, n){ int a, b, c, d; cin >> a >> b >> c >> d; L l(P(a, b), P(c, d)); rep(j, i) if(parallel(ls[j], l) && abs(cross(ls[j][0] - l[1], l[0] - l[1])) < EPS) goto SKIP; ls.pb(l); SKIP:; } n = ls.size(); if(n < 3){ cout << "Many" << endl; continue; } P p(INF, INF); G cand; bool pa[100][100] = {}; rep(i, n) rep(j, i) pa[i][j] = pa[j][i] = parallel(ls[i], ls[j]); rep(i, n) rep(j, i) if(pa[i][j]) rep(k, j) if(pa[j][k]){ cout << "None" << endl; goto END; } rep(i, n) rep(j, i) rep(k, j){ if(pa[i][j] || pa[j][k] || pa[k][i]) continue; P x = crosspoint(ls[i], ls[j]), y = crosspoint(ls[j], ls[k]), z = crosspoint(ls[k], ls[i]); if(abs(x - y) < EPS && abs(y - z) < EPS){ cand.pb(x); goto OUT; } cand.pb(incenter(x, y, z)); G g; g.pb(x); g.pb(y); g.pb(z); g = excenter(g); cand.insert(cand.end(), all(g)); goto OUT; } OUT: if(cand.size()){ rep(it, cand.size()){ bool ok = 1; double dist = distanceLP(ls[0], cand[it]); rep(i, n) if(abs(dist - distanceLP(ls[i], cand[it])) > EPS){ ok = 0; break; } if(ok){ if(p.real() != INF){ cout << "Many" << endl; goto END; } p = cand[it]; } } if(p.real() == INF){ cout << "None" << endl; goto END; } } else{ sort(all(ls)); do{ if(parallel(ls[0], ls[1]) && parallel(ls[2], ls[3])) break; }while(next_permutation(all(ls))); vector<L> t; rep(i, 2){ P m = (ls[i * 2][0] + ls[i * 2 + 1][0]) * 0.5; t.pb(L(m, m + ls[i * 2][1] - ls[i * 2][0])); } p = crosspoint(t[0], t[1]); double dist = distanceLP(ls[0], p); rep(i, 4) if(abs(distanceLP(ls[i], p) - dist) > EPS){ cout << "None" << endl; goto END; } } printf("%.9f %.9f\n", p.real(), p.imag()); END:; } return 0; }