#include<bits/stdc++.h>

 using namespace std;

 typedef long double ld;

 const int maxn=1E5+;

 const ld eps=1E-;

 const ld pi=acos(-);

 const ld inf=1E12;

 int n,m;

 ld X,Y;

 struct pt

 {

     ld x,y;

     pt(ld a=,ld b=):x(a),y(b){}

     pt operator+(const pt&A){return pt(x+A.x,y+A.y);}

     pt operator-(const pt&A){return pt(x-A.x,y-A.y);}

     pt operator*(const ld&d){return pt(x*d,y*d);}

     pt operator/(const ld&d){return pt(x/d,y/d);}

     ld operator*(const pt&A){return x*A.y-y*A.x;}

     inline void out()

     {

         cout<<"("<<x<<","<<y<<")";

     }

 }O;

 struct line

 {

     pt A,B;

     line(pt a=pt(),pt b=pt()):A(a),B(b){}

 }a[maxn];

 struct BIT

 {

     int tot;

     int t[maxn];

     inline void clear(){memset(t,,sizeof(t));}

     inline int lowbit(int x){return x&(-x);}

     inline int ask(int x){int sum=;while(x){sum+=t[x];x-=lowbit(x);}return sum;}

     inline void add(int x,int y){while(x<=tot){t[x]+=y;x+=lowbit(x);}}

 }B;

 int size;

 struct info

 {

     int type,num;

     ld x;

     info(int a=,int b=,ld c=):type(a),num(b),x(c){}

     bool operator<(const info&A)const

     {

         return x<A.x;

     }

 }tmp[maxn];

 inline pt intersection(line a,line b)

 {

     pt A=b.B-b.A,B=a.B-a.A,C=b.A-a.A;

     if(abs(A*B)<=eps)

         return pt(inf,inf);

     ld d=-(B*C)/(B*A);

     return b.A+A*d;

 }

 inline pt T(pt A)

 {

     swap(A.x,A.y);

     A.y=-A.y;

     return A;

 }

 inline pt perpendicular(pt A,line x)

 {

     pt B=T(x.B-x.A)+A;

     return intersection(line(A,B),x);

 }

 inline ld s(ld x)

 {

     return x*x;

 }

 int rk[][maxn];

 inline bool check(ld r)

 {

     size=;

     for(int i=;i<=n;++i)

     {

         pt A=perpendicular(O,a[i]);

         ld x=s(A.x-O.x)+s(A.y-O.y);

         if(x>=r*r+eps)

             continue;

         pt d=T((O-A)*sqrt(r*r/x-));

         pt P1=A+d,P2=A-d;

         ld x1=atan2(P1.y-O.y,P1.x-O.x);

         ld x2=atan2(P2.y-O.y,P2.x-O.x);

         if(!(x1<)&&!(x1>=))

         {

             x1=atan2(a[i].A.y-a[i].B.y,a[i].A.x-a[i].B.x);

             if(x1<)

                 x2=x1+pi;

             else

                 x2=x1-pi;

         }

         if(x1>x2)

             swap(x1,x2);

         tmp[++size]=info(,i,x1);

         tmp[++size]=info(,i,x2);

     }

     sort(tmp+,tmp+size+);

     for(int i=;i<=size;++i)

         rk[tmp[i].type][tmp[i].num]=i;

     B.clear();

     B.tot=size;

     long long sum=;

     for(int i=;i<=size;++i)

     {

         if(tmp[i].type==)

             B.add(i,);

         else

         {

             sum+=B.ask(i-)-B.ask(rk[][tmp[i].num]);

             B.add(rk[][tmp[i].num],-);

         }

     }

     return sum>=m;

 }

 int fa[maxn];

 int find(int x)

 {

     return x==fa[x]?x:fa[x]=find(fa[x]);

 }

 bool ok[maxn];

 inline ld dis(pt A,pt B)

 {

     return sqrt(s(A.x-B.x)+s(A.y-B.y));

 }

 inline ld get(ld r)

 {

     size=;

     for(int i=;i<=n;++i)

     {

         pt A=perpendicular(O,a[i]);

         ld x=s(A.x-O.x)+s(A.y-O.y);

         if(x>=r*r+eps)

             continue;

         ok[i]=;

         pt d=T((O-A)*sqrt(r*r/x-));

         pt P1=A+d,P2=A-d;

         ld x1=atan2(P1.y-O.y,P1.x-O.x);

         ld x2=atan2(P2.y-O.y,P2.x-O.x);

         if(abs(O.x-A.x)<=eps&&abs(O.y-A.y)<=eps)

         {

             x1=atan2(a[i].A.y-a[i].B.y,a[i].A.x-a[i].B.x);

             if(x1<)

                 x2=x1+pi;

             else

                 x2=x1-pi;

         }

         if(x1>x2)

             swap(x1,x2);

         tmp[++size]=info(,i,x1);

         tmp[++size]=info(,i,x2);

     }

     sort(tmp+,tmp+size+);

     for(int i=;i<=n;++i)

     for(int i=;i<=size;++i)

         fa[i]=rk[tmp[i].type][tmp[i].num]=i;

     fa[size+]=size+;

     ld sum=;

     for(int i=;i<=size;++i)

         if(tmp[i].type)

         {

             int l=rk[][tmp[i].num];

             fa[l]=l+,fa[i]=i+;

             while((l=find(l))<i)

             {

                 pt A=intersection(a[tmp[l].num],a[tmp[i].num]);

                 sum+=dis(intersection(a[tmp[l].num],a[tmp[i].num]),O);

                 ++l;

             }

         }

     return sum;

 }

 inline int R()

 {

     return rand()%-;

 }

 int main()

 {

 //    freopen("stigmata.in","r",stdin);

 //    freopen("stigmata.out","w",stdout);

     ios::sync_with_stdio(false);

     cin>>n>>X>>Y>>m;

     X/=,Y/=;

     for(int i=;i<=n;++i)

     {

         ld x,y;

         cin>>x>>y;

         x/=,y/=;

         a[i]=line(pt(,y),pt(,*x+y));

     }

     O=pt(X,Y);

     ld L=,R=,mid;

     while(abs(L-R)>0.00000001)

     {

         mid=(L+R)/;

         if(check(mid))

             R=mid;

         else

             L=mid;

     }

     cout<<fixed<<setprecision()<<get(R)<<endl;

     return ;

 }