

function y = callBS(S,K,T, sigma, r)
d_plus = ((r + 1/2 * sigma^2) * T - log(K/S))/ (sigma * sqrt(T) );
d_minus = ((r - 1/2 * sigma^2) * T - log(K/S))/ (sigma * sqrt(T) );

y = S * normcdf(d_plus) - K * exp(-r*T) * normcdf(d_minus);
end


