 %infs.m
for i=1:50
    s(i)=i*0.01;
end ls=length(s); m=200000;c=250;a=0.001;

Fa=zeros(1,ls); for i=1:ls
    Fa(i)=zgamma(20,m,s(i)*m,c,a)*(1+a)-c/a;
end

T=zeros(1,ls); for i=1:ls
    g(i)=zgamma(20,m,s(i)*m,c,a);
    p(i)=1/2*erfc((g(i)-m-c/a)/sqrt(2*(s(i)*m)^2));
    T(i)=(1-p(i))/p(i);
end
