题目大意：
有n条直线，求它们两两之间的交点到(p,q)前m近的距离和。
(n le 50000, m le 10^7)
题解：
二分答案r，肯定的。
接着就是求有多少个交点到(p,q)的距离<=r
每条直线可以转换成一条线段，相当于求线段交点个数。
我所知道的最优做法扫描线+平衡树+堆，复杂度：(O((m+n)logn))，其中m为交点个数。
这个这么做显然是TlE的。
注意到这些线段全部都是以(p,q)为圆心，r为半径的圆的弦。
两条弦相交相当于对应的弧真相交（相交且不包含），这里的弧不管是优弧还是劣弧都可以。
解一元二次方程求交点，atan2求出几角区间，离散化一下，就是线段树能做的事了。
注意到这题要求的是距离和。
先二分出求出最大的r使得count(r)<=m。
再考虑求出r以内的所有交点，和(p,q)求距离加入答案，最后加上(m-count(r))*r即可。
求出所有交点：
考虑再线段树的每个点上多记一个vector，表示该区间内的直线有哪些，查询时搞一下就好了。
时间复杂度：(O(n ~ log ~n ~ log ~r + (m + n ~ log ~n)))
Code:
#include<bits/stdc++.h>
#define fo(i, x, y) for(int i = x, _b = y; i <= _b; i ++)
#define ff(i, x, y) for(int i = x, _b = y; i < _b; i ++)
#define fd(i, x, y) for(int i = x, _b = y; i >= _b; i --)
#define ll long long
#define pp printf
#define hh pp("
")
using namespace std;
#define db double
const db eps = 1e-8;
const db pi = acos(-1);
const int N = 50005;
int n, m; db p, q;
db a[N], b[N];
db sqr(db x) { return x * x;}
db r;
db dis(db x, db y, db u, db v) {
return sqrt(sqr(x - u) + sqr(y - v));
}
struct P {
db x; int i;
P(db _x = 0, int _i = 0) { x = _x, i = _i;}
};
P d[N * 2]; int d0;
int cmpd(P a, P b) {
if(a.x == b.x) return a.i > b.i;
return a.x < b.x;
}
int u[N], v[N];
void build() {
d0 = 0;
fo(i, 1, n) {
db a0 = a[i] * a[i] + 1;
db b0 = -2 * p + 2 * a[i] * (b[i] - q);
db c0 = sqr(p) + sqr(b[i] - q) - sqr(r);
db t = sqr(b0) - 4 * a0 * c0;
if(t < 0) continue;
t = sqrt(t);
db x0 = (-b0 + t) / 2 / a0, x1 = (-b0 - t) / 2 / a0;
db u = atan2(a[i] * x0 + b[i] - q, x0 - p);
db v = atan2(a[i] * x1 + b[i] - q, x1 - p);
if(u > v) swap(u, v);
d[++ d0] = P(u, i);
d[++ d0] = P(v, -i);
}
sort(d + 1, d + d0 + 1, cmpd);
fo(i, 1, d0) if(d[i].i > 0)
u[d[i].i] = i; else v[-d[i].i] = i;
}
int kq;
db sum;
void solve(int x, int y) {
db x0 = (b[y] - b[x]) / (a[x] - a[y]);
sum += dis(x0, a[x] * x0 + b[x], p, q);
}
#define i0 i + i
#define i1 i + i + 1
int t[N * 8], pl, pr, px, py;
vector<int> g[N * 8];
void add(int i, int x, int y) {
if(kq) g[i].push_back(py);
t[i] ++;
if(x == y) return;
int m = x + y >> 1;
if(pl <= m) add(i0, x, m); else add(i1, m + 1, y);
t[i] = t[i0] + t[i1];
}
void ft(int i, int x, int y) {
if(y < pl || x > pr) return;
if(x >= pl && y <= pr) {
px += t[i];
if(kq) {
ff(j, 0, g[i].size())
solve(py, g[i][j]);
}
return;
}
int m = x + y >> 1;
ft(i0, x, m); ft(i1, m + 1, y);
}
ll count(db _r) {
r = _r;
build();
ll ans = 0;
fo(i, 1, d0 * 4) t[i] = 0;
fo(i, 1, d0) if(d[i].i > 0) {
pl = u[d[i].i], pr = v[d[i].i]; px = 0; py = d[i].i;
ft(1, 1, d0);
ans += px;
pl = pr = v[d[i].i];
add(1, 1, d0);
}
return ans;
}
int main() {
freopen("stigmata.in", "r", stdin);
freopen("stigmata.out", "w", stdout);
scanf("%d %lf %lf %d", &n, &p, &q, &m);
p /= 1000, q /= 1000;
fo(i, 1, n) {
scanf("%lf %lf", &a[i], &b[i]);
a[i] /= 1000, b[i] /= 1000;
}
db ans = 0;
for(db v = 1e6; v > 1e-11; v /= 2)
if(count(ans + v) <= m) ans += v;
kq = 1; int cnt = count(ans);
pp("%.9lf
", sum + ans * (m - cnt));
}