diff --git a/nifty5/__init__.py b/nifty5/__init__.py index 52d8311475f1f837f9951715b3f164585a9d62bb..c96769eccbb7411e2dd05e49f68326dbc51c22a3 100644 --- a/nifty5/__init__.py +++ b/nifty5/__init__.py @@ -82,7 +82,8 @@ from .library.dynamic_operator import (dynamic_operator, from .library.light_cone_operator import LightConeOperator from .library.wiener_filter_curvature import WienerFilterCurvature -from .library.correlated_fields import CorrelatedField, MfCorrelatedField +from .library.correlated_fields import (CorrelatedField, MfCorrelatedField, + MfPartiallyCorrelatedField) from .library.adjust_variances import (make_adjust_variances_hamiltonian, do_adjust_variances) diff --git a/nifty5/library/correlated_fields.py b/nifty5/library/correlated_fields.py index 171eba8692431e700858beca13d48d472e73ad54..a45f020d577ed1ae92f95a613334e0345e5ada84 100644 --- a/nifty5/library/correlated_fields.py +++ b/nifty5/library/correlated_fields.py @@ -141,3 +141,54 @@ def MfCorrelatedField(target, amplitudes, name='xi'): vol = reduce(mul, [sp.scalar_dvol**-0.5 for sp in hsp]) return ht(vol*A*ducktape(hsp, None, name)) + + +def MfPartiallyCorrelatedField(target, energy_amplitude, name='xi_u'): + """Constructs an operator which turns white Gaussian excitation fields + into a correlated field defined on a DomainTuple with two entries. + On the first domain, a white correlation structure is assumed. On + the second domain, the correlation structures given by energy_amplitude + is used. + + This operator may be used as a model for multi-frequency reconstructions + with correlation structure only in the energy direction. + + Parameters + ---------- + target : Domain, DomainTuple or tuple of Domain + Target of the operator. Must contain exactly two spaces. + It is assumed that the second space is the energy domain. + energy_amplitude: Operator + amplitude operator for the energy correlation structure + name : string + :class:`MultiField` key for xi-field. + + Returns + ------- + Operator + Correlated field + + Notes + ----- + In NIFTy, non-harmonic RGSpaces are by definition periodic. Therefore + the operator constructed by this method will output a correlated field + with *periodic* boundary conditions. If a non-periodic field is needed, + one needs to combine this operator with a :class:`FieldZeroPadder` or even + two (one for the energy and one for the spatial subdomain) + """ + tgt = DomainTuple.make(target) + if len(tgt) != 2: + raise ValueError + + h_space = DomainTuple.make([dom.get_default_codomain() for dom in tgt]) + ht1 = HarmonicTransformOperator(h_space, target=tgt[0], space=0) + ht2 = HarmonicTransformOperator(ht1.target, target=tgt[1], space=1) + ht = ht2 @ ht1 + + p_space = energy_amplitude.target[0] + power_distributor = PowerDistributor(h_space[-1], p_space) + A = power_distributor(energy_amplitude) + + dd = ContractionOperator(h_space, 0).adjoint + + return ht((dd @ A)*ducktape(h_space, None, name))