Add Uniform distribution operator

nifty5/__init__.py

@@ -76,7 +76,8 @@ from .plot import Plot
 
 from .library.smooth_linear_amplitude import (
     SLAmplitude, LinearSLAmplitude, CepstrumOperator)
-from .library.inverse_gamma_operator import InverseGammaOperator
+from .library.probability_transform_operators import (InverseGammaOperator,
+                                                      UniformOperator)
 from .library.los_response import LOSResponse
 from .library.dynamic_operator import (dynamic_operator,
                                        dynamic_lightcone_operator)

nifty5/library/inverse_gamma_operator.py → nifty5/library/probability_transform_operators.py

@@ -104,3 +104,45 @@ class InverseGammaOperator(Operator):
     @property
     def q(self):
         return self._q
+
+
+class UniformOperator(Operator):
+    """Transforms a Gaussian into a uniform distribution in the interval [0,1].
+
+    Parameters
+    ----------
+    domain : Domain, tuple of Domain or DomainTuple
+        The domain on which the field shall be defined. This is at the same
+        time the domain and the target of the operator.
+    scale : float
+        Scales the output, default: 1.
+    offset : float
+        Shifts the output. The shift is applied after the scaling, default: 0.
+    """
+    def __init__(self, domain, scale=1, offset=0):
+        self._domain = self._target = DomainTuple.make(domain)
+        self._scale, self._offset = float(scale), float(offset)
+        if self._scale == 0.:
+            raise ValueError
+
+    def apply(self, x):
+        self._check_input(x)
+        lin = isinstance(x, Linearization)
+
+        val = x.val.local_data if lin else x.local_data
+        res = norm.cdf(val)
+        if self._scale != 1.:
+            res *= self._scale
+        if self._offset != 0.:
+            res += self._offset
+
+        resfld = Field.from_local_data(self._domain, res)
+        if not lin:
+            return resfld
+
+        der = norm.pdf(val)
+        if self._scale != 1.:
+            der *= self._scale
+        jac = makeOp(Field.from_local_data(self._domain, der))
+        jac = jac(x.jac)
+        return x.new(resfld, jac)

test/test_operators/test_jacobian.py

@@ -130,6 +130,15 @@ def testPointModel(space, seed):
     ift.extra.check_jacobian_consistency(model, pos, tol=1e-2, ntries=20)
 
 
+@pmp('scale', [1, 2.2])
+@pmp('offset', [0, 127.46])
+def testUniformOperator(space, seed, scale, offset):
+    S = ift.ScalingOperator(1., space)
+    pos = S.draw_sample()
+    model = ift.UniformOperator(space, scale, offset)
+    ift.extra.check_jacobian_consistency(model, pos, ntries=20)
+
+
 @pmp('target', [
     ift.RGSpace(64, distances=.789, harmonic=True),
     ift.RGSpace([32, 32], distances=.789, harmonic=True),