sedFitThin | R Documentation |
Function takes Herschel-SPIRE photometry and fits optically-thin greybody function for a single-component temperature and galaxy luminosity. Function generates nsamp realizations of observed flux densities with standard deviations for error analysis.
sedFitThin(s, e = s*0.2, z = 2.5, nsamp = 100, alpha = 2, beta = 1.5, wl= c(250, 350, 500), sc.start = 1.e-6, T.start = 50)
s |
Vector of observed-frame flux densities [Jy] |
e |
Vector of standard deviation of observed-frame flux density [Jy] |
z |
Galaxy redshift |
nsamp |
Number of realizations for Monte-Carlo calculation |
alpha |
Index of power-law for short-wavelength extension |
beta |
Dust emissivity power law |
wl |
Vector of observed-frame wavelengths corresponding to |
sc.start |
Initial guess for fit luminosity density scaling factor |
T.start |
Initial guess for dust temperature [K] |
Conversion from observed to rest frame is from equation (24) in Hogg 2000. Dust temperature and 8-1000 micron luminosity derivation is described in Blain, Barnard & Chapman 2003. Galaxy SEDs typically fall off more slowly than greybody on the Wien side; see plot generated by examples below to visualize power-law extension suggested by Blain et al. 2003.
List of class sedfit
with elements:
td |
Mean of dust temperature distribution |
e.td |
Standard deviation of dust temperature distribution |
lum.gb |
Mean of greybody luminosity distribution |
e.lum.gb |
Standard deviation of greybody luminosity distribution |
lum.gbpl |
Mean of greybody-power law luminosity distribution |
e.lum.gbpl |
Standard deviation of greybody-power law luminosity distribution |
scaleFactor |
Conversion between observed frame flux density and rest frame luminosity density |
success |
Fraction of fit attempts that converged |
results |
Matrix with |
Fit will sometimes crash on numerical derivative and throw an error. In
this case the routine will halt without producing results. The more
usual lack of convergence is reported as a warning, and the
corresponding results will be NA
in the output matrix.
A. Harris
Hogg 2000, astro-ph 9905116v4; Blain, Barnard & Chapman 2003, MNRAS 338, 733.
s <- c(0.242, 0.293, 0.231) e <- c(0.037, 0.045, 0.036) z <- 2.41 beta <- 1.5 alpha <- 2 X <- sedFitThin(s=s, e=e, z=z, alpha=alpha, beta=beta, nsamp=100) str(X) ## Make a plot # greybody in blue, power-law extension in red dashed line # functions # optically thin greybody otGreybody <- function(nu, T, beta, sc=1) { # nu in GHz, T in K, beta and sc unitless sc*nu^(3+beta)/(exp(0.04801449*nu/T) - 1) } # high frequency tail hfTail <- function(nu, alpha) nu^-alpha # # setups for 8-1000 microns: nu.low <- 3e5/1000 nu.high <- 3e5/8 l.nue <- s*X$scaleFactor # # greybody nue.sweep <- seq(nu.low, nu.high, len=350) pred <- otGreybody(nue.sweep, X$results[1,1], beta=beta, X$results[1,4]) ylim <- range(pred, l.nue) par(fig=c(0,1,0.2,1), mgp=c(1.8, 0.6, 0)) plot(3e5/nue.sweep, pred, t='l', ylim=ylim, log='xy', col=4, xlab='Rest frame wavelength [microns]', ylab=expression(paste('Luminosity density [ ', L[sun], ' ', Hz^-1, ']'))) # power law nue.sweep <- seq(X$results[1,5], nu.high, len=100) val.t <- otGreybody(nu=X$results[1,5], T=X$results[1,1], beta=beta, sc=X$results[1,4]) lines(3e5/nue.sweep, val.t*hfTail(nue.sweep/X$results[1,5], alpha=alpha), col=2, lwd=1, lty=2) # data wl <- c(250, 350, 500) nue <- 3e5/wl*(1+z) points(3e5/nue, l.nue, pch=16, col=3)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.