def lnlike_spec(spec_mu, obs=None, spec_noise=None, **vectors):
"""Calculate the likelihood of the spectroscopic data given the
spectroscopic model. Allows for the use of a gaussian process
covariance matrix for multiplicative residuals.
:param spec_mu:
The mean model spectrum, in linear or logarithmic units, including
e.g. calibration and sky emission.
:param obs: (optional)
A dictionary of the observational data, including the keys
*``spectrum`` a numpy array of the observed spectrum, in linear or
logarithmic units (same as ``spec_mu``).
*``unc`` the uncertainty of same length as ``spectrum``
*``mask`` optional boolean array of same length as ``spectrum``
*``wavelength`` if using a GP, the metric that is used in the
kernel generation, of same length as ``spectrum`` and typically
giving the wavelength array.
:param spec_noise: (optional)
A NoiseModel object with the methods `compute` and `lnlikelihood`.
If ``spec_noise`` is supplied, the `wavelength` entry in the obs
dictionary must exist.
:param vectors: (optional)
A dictionary of vectors of same length as ``wavelength`` giving
possible weghting functions for the kernels
:returns lnlikelhood:
The natural logarithm of the likelihood of the data given the mean
model spectrum.
"""
if obs['spectrum'] is None:
return 0.0
mask = obs.get('mask', slice(None))
vectors['mask'] = mask
vectors['wavelength'] = obs['wavelength']
delta = (obs['spectrum'] - spec_mu)[mask]
if spec_noise is not None:
try:
spec_noise.compute(**vectors)
return spec_noise.lnlikelihood(delta)
except(LinAlgError):
return np.nan_to_num(-np.inf)
else:
# simple noise model
var = (obs['unc'][mask])**2
lnp = -0.5*( (delta**2/var).sum() + np.log(2*np.pi*var).sum() )
return lnp
评论列表
文章目录
