function y0 = ft(n,pt)

if pt < 0 || pt > 1
    disp('输入有误，请输入[0 1]区间的数值')
elseif 0.5 <= pt && pt < 1
    p1 = 2 * (1-pt);
    tp = sqrt(n/fbeta(n,p1) - n);
elseif 0 < pt && pt < 0.5
    p2 = 2*pt;
    tp = -sqrt(n/fbeta(n,p2) - n);
elseif pt == 1
    tp = inf;
else
    tp = 0;
end
y0 = tp;

%beta分布的分位数计算
function y1=fbeta(n,p)

b = 1;
a = 0;
if abs(f(b,n,p)) < 1e-10
    y1 = b;
else
    while abs(b - a) > 1e-6
     tb = (b + a)/2;
     if abs(f(tb,n,p)) < 1e-10
         y1 = tb;
         break;
     elseif f(tb,n,p) * f(b,n,p) < 0
        a = tb;
     else
        b = tb;        
     end     
    end
    y1 = tb;
end


%beta分布方程f = I(a,b) - p
function y2=f(z,n,p)

syms x t;
f1=int(t^(n/2-1)*(1-t)^(-0.5),0,'x');
f2=f1/beta(n/2,0.5);
f3=f2-p;
f3=subs(f3,x,z);
y2=eval('f3');