R code used to produce the graph shown below.

airquality$Month <- factor(airquality$Month) 
model <- lm(Temp ~ Month - 1, data=airquality) 
mean.values <- coef(summary(model))[,1] 
mean.values 
 
sem <- coef(summary(model))[,2] 
scale <- qt(0.975, df=summary(model)$df[2]) 
barplot2(mean.values, col="blue", 
         xlab="Month", ylab="Temperature (F)", 
         plot.ci=TRUE, 
         ci.u = mean.values+scale*sem, 
         ci.l = mean.values) 
 
airquality$Month <- factor(airquality$Month) 
model <- lm(Temp ~ Month - 1, data=airquality) 
mean.values <- coef(summary(model))[,1] 
sem <- coef(summary(model))[,2] 
scale <- qt(0.975, df=summary(model)$df[2]) 
barplot2(mean.values, col=myblue,  ylim=c(0,99), 
         xlab="Month", ylab="Temperature (F)", 
         plot.ci=TRUE, 
         ci.u = mean.values+scale*sem, 
         ci.l = mean.values, ci.lwd=3)