[HNOI2019]鱼（计算几何）

看到数据范围n<=1000，但感觉用O(n^2)不现实，所以考虑方向应该是O(n^2logn)。

一种暴力做法：用vector存到1点相同的2点和到2点相同的1点，然后枚举A，枚举BC，再枚举D，然后枚举EF，O(n^4)，但复杂度远远不满（符合条件的太少，而且也卡不掉），所以可以获得40pts的好成绩。

正解：

首先可以想到一点：确定A,D两点后，因为鱼身和鱼尾互不影响，可以分开计算鱼身和鱼尾，然后将二者数量相乘即可。鱼尾比较好算，直接把所有点按照极角坐标排序，然后枚举A时转移即可，复杂度O(n^2)；鱼身挺难算的，因为BC中点要在AD上，且与AD垂直（即斜率相乘为-1，或一个平行一个竖直），这时候，就要预处理两两点对的中心和斜率，然后查找时二分找出范围。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1100;
const double pi=acos(-1),eps=1e-10;
struct point{
    ll x,y;double k;
    void cal(){k=atan2(y,x);}
    bool operator<(const point&b)const{return k<b.k;}
}a[N],q[N<<1];
struct node{
    ll a,b,c,d;
    bool operator<(const node &t)const{
        if(a!=t.a)return a<t.a;
        if(b!=t.b)return b<t.b;
        if(c!=t.c)return c<t.c;
        return d<t.d;
    }
}s[N*N];
int n,m,tot,len,num[N];
ll ans,sum,dis[N],b[N];
void add(int x){sum+=num[dis[(x-1)%m+1]],num[dis[(x-1)%m+1]]++;}
void del(int x){num[dis[(x-1)%m+1]]--,sum-=num[dis[(x-1)%m+1]];}
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++)scanf("%lld%lld",&a[i].x,&a[i].y);
    for(int i=1;i<=n;i++)
    for(int j=i+1;j<=n;j++)
    {
        ll x=a[i].x-a[j].x,y=a[i].y-a[j].y,d=abs(__gcd(x,y));
        x/=d,y/=d;
        if(x<0)x=-x,y=-y;
        if(!x)y=abs(y);
        s[++tot]=(node){x,y,x*(a[i].x+a[j].x)+y*(a[i].y+a[j].y),x*a[i].y-y*a[i].x};
    }
    sort(s+1,s+tot+1);
    for(int i=1;i<=n;i++)
    {
        m=0;
        for(int j=1;j<=n;j++)if(i!=j)q[++m]=(point){a[j].x-a[i].x,a[j].y-a[i].y},q[m].cal();
        sort(q+1,q+m+1);
        for(int j=m+1;j<=2*m;j++)q[j]=q[j-m],q[j].k+=2*pi;
        len=0;
        for(int j=1;j<=m;j++)dis[j]=b[++len]=q[j].x*q[j].x+q[j].y*q[j].y;
        sort(b+1,b+len+1);
        len=unique(b+1,b+len+1)-b-1;
        for(int j=1;j<=m;j++)dis[j]=lower_bound(b+1,b+len+1,dis[j])-b;
        sum=0;
        memset(num,0,sizeof num);
        for(int j=1,qs=0,qe=0;j<=m;j++)
        {
            while(q[qs+1].k+eps<q[j].k+1.5*pi)add(++qs);
            while(q[qe+1].k-eps<q[j].k+0.5*pi)del(++qe);
            ll x=-q[j].y,y=q[j].x,d=abs(__gcd(x,y));
            x/=d,y/=d;
            if(x<0)x=-x,y=-y;
            if(!x)y=abs(y);
            ll t=x*(q[j].x+a[i].x+a[i].x)+y*(q[j].y+a[i].y+a[i].y);
            ll L=x*a[i].y-y*a[i].x,R=x*(q[j].y+a[i].y)-y*(q[j].x+a[i].x);
            if(L>R)swap(L,R);
            node xl=(node){x,y,t,L},xr=(node){x,y,t,R};
            ans+=sum*(lower_bound(s+1,s+tot+1,xr)-upper_bound(s+1,s+tot+1,xl));
        }
    }
    cout<<ans*4;
}