diff --git a/nifty4/library/__init__.py b/nifty4/library/__init__.py index 1efb22aef60e24e20640ffd6f0d12025d91be409..621ca6119053525283ec833f5cf75cdca8473496 100644 --- a/nifty4/library/__init__.py +++ b/nifty4/library/__init__.py @@ -5,5 +5,6 @@ from .nonlinear_wiener_filter_energy import NonlinearWienerFilterEnergy from .nonlinearities import Exponential, Linear, PositiveTanh, Tanh from .poisson_energy import PoissonEnergy from .unit_log_gauss import UnitLogGauss +from .poisson_log_likelihood import PoissonLogLikelihood from .wiener_filter_curvature import WienerFilterCurvature from .wiener_filter_energy import WienerFilterEnergy diff --git a/nifty4/library/poisson_log_likelihood.py b/nifty4/library/poisson_log_likelihood.py new file mode 100644 index 0000000000000000000000000000000000000000..4ddb041d9318f2dc9a42c2201f71f25e576631fd --- /dev/null +++ b/nifty4/library/poisson_log_likelihood.py @@ -0,0 +1,44 @@ +from numpy import inf, isnan + +from ..minimization import Energy +from ..operators import SandwichOperator +from ..sugar import log, makeOp + + +class PoissonLogLikelihood(Energy): + def __init__(self, position, lamb, d): + """ + s: Sky model object + + value = 0.5 * s.vdot(s), i.e. a log-Gauss distribution with unit + covariance + """ + super(PoissonLogLikelihood, self).__init__(position) + self._lamb = lamb.at(position) + self._d = d + + lamb_val = self._lamb.value + + self._value = lamb_val.sum() - d.vdot(log(lamb_val)) + if isnan(self._value): + self._value = inf + self._gradient = self._lamb.gradient.adjoint_times(1 - d/lamb_val) + + # metric = makeOp(d/lamb_val/lamb_val) + metric = makeOp(1./lamb_val) + self._curvature = SandwichOperator.make(self._lamb.gradient, metric) + + def at(self, position): + return self.__class__(position, self._lamb, self._d) + + @property + def value(self): + return self._value + + @property + def gradient(self): + return self._gradient + + @property + def curvature(self): + return self._curvature