Hi Tim,
If you apply the usual large sample symmetric interval then sure, you could obtain CIs. The reason that this isn't done by default is that many factor score computations include some Bayesian prior information (often borrowed from the model itself), so you don't really get CIs but rather the Bayesian analogue symmetric "credible regions", in which case the usual large-sample quantiles are not really kosher anymore. See below on how to compute CIs for ML estimates, in which this strategy is more in-line with theory.
#---------------------------
library(mirt)
mod <- mirt(Science, 1)
head(fs)
# 95% large-sample symmetric CIs: x +- q * SE
CIs <- t(apply(fs, 1, function(x) x['F1'] + qnorm(c(.025, .975)) * x['SE_F1']))
head(CIs)
#---------------------------