function y=BS_call(S0,K,T,r,sigma)
    tic();
    d1=(log(S0/K)+(r+sigma^2/2)*T)/(sigma*sqrt(T));
    d2=(log(S0/K)+(r-sigma^2/2)*T)/(sigma*sqrt(T));
    price=S0*cdfnor("PQ",d1,0,1)-K*exp(-r*T)*cdfnor("PQ",d2,0,1);    
    time=toc();
    y=[price time]
endfunction