function [d,g]=digitaldelta(so,x,r,sig) 
% like blsdelta only for a digital option which pays $1  if S(1)>x  and
% otherwise pays $0. Outputs both the delta d and the gamma g.
%corrected Jan 2005
d1 = (log(so./x)+(r-sig.^2/2))./sig;
d=normpdf(d1)/(x*sig);
g=-d.*(d1./(so*sig));