You need to use ?segments. Try this:
x.vals <- c(mean(g)-(1:3*sd(g)), mean(g), mean(g)+(1:3*sd(g)))
y.vals <- dnorm(x.vals, mean=mean(g),sd=sd(g))*diff(h$mids[1:2])*length(g)
segments(x0=x.vals, y0=0, x1=x.vals, y1=y.vals)
Note that not all SDs show up on the plot. To make the original plot wider, you could use xlim=c(-3,13)
in the call to hist()
(and you'd also need to change the range of xfit
accordingly).