field.py 22.7 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
Theo Steininger's avatar
Theo Steininger committed
13
14
15
16
17
#
# Copyright(C) 2013-2017 Max-Planck-Society
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
# and financially supported by the Studienstiftung des deutschen Volkes.
18

Martin Reinecke's avatar
Martin Reinecke committed
19
from __future__ import division, print_function
Martin Reinecke's avatar
Martin Reinecke committed
20
from builtins import range
csongor's avatar
csongor committed
21
import numpy as np
Martin Reinecke's avatar
Martin Reinecke committed
22
23
24
from .spaces.power_space import PowerSpace
from . import nifty_utilities as utilities
from .random import Random
Martin Reinecke's avatar
Martin Reinecke committed
25
from .domain_tuple import DomainTuple
Martin Reinecke's avatar
Martin Reinecke committed
26
from functools import reduce
27
from . import dobj
28

csongor's avatar
csongor committed
29

Martin Reinecke's avatar
Martin Reinecke committed
30
class Field(object):
Theo Steininger's avatar
Theo Steininger committed
31
32
33
    """ The discrete representation of a continuous field over multiple spaces.

    In NIFTY, Fields are used to store data arrays and carry all the needed
34
    metainformation (i.e. the domain) for operators to be able to work on them.
Martin Reinecke's avatar
updates    
Martin Reinecke committed
35
    In addition, Field has methods to work with power spectra.
Theo Steininger's avatar
Theo Steininger committed
36

37
38
39
40
    Parameters
    ----------
    domain : DomainObject
        One of the space types NIFTY supports. RGSpace, GLSpace, HPSpace,
Theo Steininger's avatar
Theo Steininger committed
41
        LMSpace or PowerSpace. It might also be a FieldArray, which is
42
        an unstructured domain.
Theo Steininger's avatar
Theo Steininger committed
43

Martin Reinecke's avatar
stage1    
Martin Reinecke committed
44
    val : scalar, numpy.ndarray, Field
45
46
47
        The values the array should contain after init. A scalar input will
        fill the whole array with this scalar. If an array is provided the
        array's dimensions must match the domain's.
Theo Steininger's avatar
Theo Steininger committed
48

49
    dtype : type
Martin Reinecke's avatar
updates    
Martin Reinecke committed
50
        A numpy.type. Most common are float and complex.
Theo Steininger's avatar
Theo Steininger committed
51

52
53
54
55
    copy: boolean

    Attributes
    ----------
Martin Reinecke's avatar
stage1    
Martin Reinecke committed
56
    val : numpy.ndarray
Theo Steininger's avatar
Theo Steininger committed
57

Martin Reinecke's avatar
Martin Reinecke committed
58
    domain : DomainTuple
59
60
61
        See Parameters.
    dtype : type
        Contains the datatype stored in the Field.
Theo Steininger's avatar
Theo Steininger committed
62

63
64
65
66
67
68
69
    Raise
    -----
    TypeError
        Raised if
            *the given domain contains something that is not a DomainObject
             instance
            *val is an array that has a different dimension than the domain
Theo Steininger's avatar
Theo Steininger committed
70

71
    """
72

Theo Steininger's avatar
Theo Steininger committed
73
    # ---Initialization methods---
74

Martin Reinecke's avatar
stage1    
Martin Reinecke committed
75
    def __init__(self, domain=None, val=None, dtype=None, copy=False):
76
        self.domain = self._parse_domain(domain=domain, val=val)
77

Martin Reinecke's avatar
Martin Reinecke committed
78
        dtype = self._infer_dtype(dtype=dtype, val=val)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
79
        if isinstance(val, Field):
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
80
            if self.domain != val.domain:
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
81
                raise ValueError("Domain mismatch")
82
            self._val = dobj.from_object(val.val, dtype=dtype, copy=copy)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
83
        elif (np.isscalar(val)):
84
85
            self._val = dobj.full(self.domain.shape, dtype=dtype, fill_value=val)
        elif isinstance(val, dobj.data_object):
Martin Reinecke's avatar
Martin Reinecke committed
86
            if self.domain.shape == val.shape:
87
                self._val = dobj.from_object(val, dtype=dtype, copy=copy)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
88
89
90
            else:
                raise ValueError("Shape mismatch")
        elif val is None:
91
            self._val = dobj.empty(self.domain.shape, dtype=dtype)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
92
93
        else:
            raise TypeError("unknown source type")
csongor's avatar
csongor committed
94

Martin Reinecke's avatar
Martin Reinecke committed
95
96
    @staticmethod
    def _parse_domain(domain, val=None):
97
        if domain is None:
98
            if isinstance(val, Field):
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
99
100
                return val.domain
            if np.isscalar(val):
101
                return DomainTuple.make(())  # empty domain tuple
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
102
            raise TypeError("could not infer domain from value")
Martin Reinecke's avatar
Martin Reinecke committed
103
        return DomainTuple.make(domain)
104

Martin Reinecke's avatar
Martin Reinecke committed
105
    # MR: this needs some rethinking ... do we need to have at least float64?
Martin Reinecke's avatar
Martin Reinecke committed
106
107
    @staticmethod
    def _infer_dtype(dtype, val):
Martin Reinecke's avatar
Martin Reinecke committed
108
109
110
111
112
        if val is None or dtype is not None:
            return np.result_type(dtype, np.float64)
        if isinstance(val, Field):
            return val.dtype
        return np.result_type(val, np.float64)
113

114
    # ---Factory methods---
115

Martin Reinecke's avatar
Martin Reinecke committed
116
117
    @staticmethod
    def from_random(random_type, domain, dtype=np.float64, **kwargs):
118
119
120
121
122
123
124
        """ Draws a random field with the given parameters.

        Parameters
        ----------
        random_type : String
            'pm1', 'normal', 'uniform' are the supported arguments for this
            method.
Theo Steininger's avatar
Theo Steininger committed
125

126
127
        domain : DomainObject
            The domain of the output random field
Theo Steininger's avatar
Theo Steininger committed
128

129
130
        dtype : type
            The datatype of the output random field
Theo Steininger's avatar
Theo Steininger committed
131

132
133
134
135
136
137
138
        Returns
        -------
        out : Field
            The output object.

        See Also
        --------
139
        power_synthesize
140
        """
Theo Steininger's avatar
Theo Steininger committed
141

Martin Reinecke's avatar
Martin Reinecke committed
142
        domain = DomainTuple.make(domain)
143
        generator_function = getattr(Random, random_type)
144
        return Field(domain=domain, val=generator_function(dtype=dtype,
Martin Reinecke's avatar
Martin Reinecke committed
145
                     shape=domain.shape, **kwargs))
146

147
148
    # ---Powerspectral methods---

Martin Reinecke's avatar
Martin Reinecke committed
149
150
    def power_analyze(self, spaces=None, binbounds=None,
                      keep_phase_information=False):
Theo Steininger's avatar
Theo Steininger committed
151
        """ Computes the square root power spectrum for a subspace of `self`.
Theo Steininger's avatar
Theo Steininger committed
152

Theo Steininger's avatar
Theo Steininger committed
153
154
155
        Creates a PowerSpace for the space addressed by `spaces` with the given
        binning and computes the power spectrum as a Field over this
        PowerSpace. This can only be done if the subspace to  be analyzed is a
156
        harmonic space. The resulting field has the same units as the initial
Theo Steininger's avatar
Theo Steininger committed
157
        field, corresponding to the square root of the power spectrum.
158
159
160

        Parameters
        ----------
Theo Steininger's avatar
Theo Steininger committed
161
        spaces : int *optional*
Martin Reinecke's avatar
Martin Reinecke committed
162
            The subspace for which the powerspectrum shall be computed.
Theo Steininger's avatar
Theo Steininger committed
163
164
165
            (default : None).
        binbounds : array-like *optional*
            Inner bounds of the bins (default : None).
Martin Reinecke's avatar
Martin Reinecke committed
166
            if binbounds==None : bins are inferred.
167
168
169
170
171
172
173
174
175
176
        keep_phase_information : boolean, *optional*
            If False, return a real-valued result containing the power spectrum
            of the input Field.
            If True, return a complex-valued result whose real component
            contains the power spectrum computed from the real part of the
            input Field, and whose imaginary component contains the power
            spectrum computed from the imaginary part of the input Field.
            The absolute value of this result should be identical to the output
            of power_analyze with keep_phase_information=False.
            (default : False).
Theo Steininger's avatar
Theo Steininger committed
177

178
179
        Raise
        -----
Martin Reinecke's avatar
Martin Reinecke committed
180
181
        TypeError
            Raised if any of the input field's domains is not harmonic
Theo Steininger's avatar
Theo Steininger committed
182

183
184
        Returns
        -------
Theo Steininger's avatar
Theo Steininger committed
185
        out : Field
Martin Reinecke's avatar
Martin Reinecke committed
186
            The output object. Its domain is a PowerSpace and it contains
187
188
189
190
191
            the power spectrum of 'self's field.

        See Also
        --------
        power_synthesize, PowerSpace
Theo Steininger's avatar
Theo Steininger committed
192

193
        """
Theo Steininger's avatar
Theo Steininger committed
194

Theo Steininger's avatar
Theo Steininger committed
195
        # check if all spaces in `self.domain` are either harmonic or
196
197
198
        # power_space instances
        for sp in self.domain:
            if not sp.harmonic and not isinstance(sp, PowerSpace):
Martin Reinecke's avatar
Martin Reinecke committed
199
200
                print("WARNING: Field has a space in `domain` which is "
                      "neither harmonic nor a PowerSpace.")
201
202

        # check if the `spaces` input is valid
203
        if spaces is None:
Martin Reinecke's avatar
Martin Reinecke committed
204
            spaces = range(len(self.domain))
Martin Reinecke's avatar
Martin Reinecke committed
205
206
        else:
            spaces = utilities.cast_iseq_to_tuple(spaces)
207
208

        if len(spaces) == 0:
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
209
            raise ValueError("No space for analysis specified.")
210

211
        if keep_phase_information:
212
            parts = [self.real*self.real, self.imag*self.imag]
213
        else:
214
            parts = [self.real*self.real + self.imag*self.imag]
215
216

        for space_index in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
217
218
            parts = [self._single_power_analyze(field=part,
                                                idx=space_index,
219
                                                binbounds=binbounds)
220
                     for part in parts]
221

222
        return parts[0] + 1j*parts[1] if keep_phase_information else parts[0]
223

Martin Reinecke's avatar
Martin Reinecke committed
224
    @staticmethod
Martin Reinecke's avatar
Martin Reinecke committed
225
226
227
    def _single_power_analyze(field, idx, binbounds):
        power_domain = PowerSpace(field.domain[idx], binbounds)
        pindex = power_domain.pindex
Martin Reinecke's avatar
Martin Reinecke committed
228
        axes = field.domain.axes[idx]
Martin Reinecke's avatar
Martin Reinecke committed
229
230
231
232
        new_pindex_shape = [1] * len(field.shape)
        for i, ax in enumerate(axes):
            new_pindex_shape[ax] = pindex.shape[i]
        pindex = np.broadcast_to(pindex.reshape(new_pindex_shape), field.shape)
Theo Steininger's avatar
Theo Steininger committed
233

Martin Reinecke's avatar
Martin Reinecke committed
234
        power_spectrum = utilities.bincount_axis(pindex, weights=field.val,
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
235
                                                 axis=axes)
Martin Reinecke's avatar
Martin Reinecke committed
236
237
238
239
240
241
        new_rho_shape = [1] * len(power_spectrum.shape)
        new_rho_shape[axes[0]] = len(power_domain.rho)
        power_spectrum /= power_domain.rho.reshape(new_rho_shape)
        result_domain = list(field.domain)
        result_domain[idx] = power_domain
        return Field(result_domain, power_spectrum)
242

Martin Reinecke's avatar
Martin Reinecke committed
243
    def _compute_spec(self, spaces):
244
245
        if spaces is None:
            spaces = range(len(self.domain))
Martin Reinecke's avatar
Martin Reinecke committed
246
247
        else:
            spaces = utilities.cast_iseq_to_tuple(spaces)
248
249
250
251
252
253
254
255
256

        # create the result domain
        result_domain = list(self.domain)
        for i in spaces:
            if not isinstance(self.domain[i], PowerSpace):
                raise ValueError("A PowerSpace is needed for field "
                                 "synthetization.")
            result_domain[i] = self.domain[i].harmonic_partner

257
        spec = dobj.sqrt(self.val)
258
        for i in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
259
260
261
262
            power_space = self.domain[i]
            local_blow_up = [slice(None)]*len(spec.shape)
            # it is important to count from behind, since spec potentially
            # grows with every iteration
Martin Reinecke's avatar
Martin Reinecke committed
263
            index = self.domain.axes[i][0]-len(self.shape)
Martin Reinecke's avatar
Martin Reinecke committed
264
265
266
267
            local_blow_up[index] = power_space.pindex
            # here, the power_spectrum is distributed into the new shape
            spec = spec[local_blow_up]
        return Field(result_domain, val=spec)
268
269

    def power_synthesize(self, spaces=None, real_power=True, real_signal=True):
Theo Steininger's avatar
Theo Steininger committed
270
        """ Yields a sampled field with `self`**2 as its power spectrum.
Theo Steininger's avatar
Theo Steininger committed
271

Theo Steininger's avatar
Theo Steininger committed
272
        This method draws a Gaussian random field in the harmonic partner
Martin Reinecke's avatar
typos    
Martin Reinecke committed
273
        domain of this field's domains, using this field as power spectrum.
Theo Steininger's avatar
Theo Steininger committed
274

275
276
277
        Parameters
        ----------
        spaces : {tuple, int, None} *optional*
Theo Steininger's avatar
Theo Steininger committed
278
279
280
            Specifies the subspace containing all the PowerSpaces which
            should be converted (default : None).
            if spaces==None : Tries to convert the whole domain.
281
        real_power : boolean *optional*
Theo Steininger's avatar
Theo Steininger committed
282
283
            Determines whether the power spectrum is treated as intrinsically
            real or complex (default : True).
284
        real_signal : boolean *optional*
Theo Steininger's avatar
Theo Steininger committed
285
286
            True will result in a purely real signal-space field
            (default : True).
Theo Steininger's avatar
Theo Steininger committed
287

288
289
290
291
        Returns
        -------
        out : Field
            The output object. A random field created with the power spectrum
Theo Steininger's avatar
Theo Steininger committed
292
            stored in the `spaces` in `self`.
293

Theo Steininger's avatar
Theo Steininger committed
294
295
296
297
298
299
        Notes
        -----
        For this the spaces specified by `spaces` must be a PowerSpace.
        This expects this field to be the square root of a power spectrum, i.e.
        to have the unit of the field to be sampled.

300
301
302
        See Also
        --------
        power_analyze
Theo Steininger's avatar
Theo Steininger committed
303
304
305
306
307

        Raises
        ------
        ValueError : If domain specified by `spaces` is not a PowerSpace.

308
        """
Theo Steininger's avatar
Theo Steininger committed
309

Martin Reinecke's avatar
Martin Reinecke committed
310
        spec = self._compute_spec(spaces)
311
312
313

        # create random samples: one or two, depending on whether the
        # power spectrum is real or complex
314
        result = [self.from_random('normal', mean=0., std=1.,
Martin Reinecke's avatar
Martin Reinecke committed
315
                                   domain=spec.domain,
316
317
318
319
320
321
322
                                   dtype=np.float if real_signal
                                   else np.complex)
                  for x in range(1 if real_power else 2)]

        # MR: dummy call - will be removed soon
        if real_signal:
            self.from_random('normal', mean=0., std=1.,
Martin Reinecke's avatar
Martin Reinecke committed
323
                             domain=spec.domain, dtype=np.float)
324
325

        # apply the rescaler to the random fields
326
        result[0] *= spec.real
327
        if not real_power:
328
            result[1] *= spec.imag
329

330
        return result[0] if real_power else result[0] + 1j*result[1]
331

332
    def power_synthesize_special(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
333
        spec = self._compute_spec(spaces)
334
335
336

        # MR: dummy call - will be removed soon
        self.from_random('normal', mean=0., std=1.,
Martin Reinecke's avatar
Martin Reinecke committed
337
                         domain=spec.domain, dtype=np.complex)
338
339

        return spec.real
340

Theo Steininger's avatar
Theo Steininger committed
341
    # ---Properties---
342

Theo Steininger's avatar
Theo Steininger committed
343
344
    @property
    def val(self):
Martin Reinecke's avatar
stage1    
Martin Reinecke committed
345
        """ Returns the data object associated with this Field.
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
346
        No copy is made.
Theo Steininger's avatar
Theo Steininger committed
347

348
349
        Returns
        -------
Martin Reinecke's avatar
stage1    
Martin Reinecke committed
350
        out : numpy.ndarray
351
        """
Martin Reinecke's avatar
Martin Reinecke committed
352
        return self._val
csongor's avatar
csongor committed
353

Martin Reinecke's avatar
Martin Reinecke committed
354
355
356
357
    @property
    def dtype(self):
        return self._val.dtype

358
359
    @property
    def shape(self):
Theo Steininger's avatar
Theo Steininger committed
360
        """ Returns the total shape of the Field's data array.
Theo Steininger's avatar
Theo Steininger committed
361

362
363
364
        Returns
        -------
        out : tuple
Martin Reinecke's avatar
Martin Reinecke committed
365
            The output object. The tuple contains the dimensions of the spaces
366
            in domain.
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
367
       """
Martin Reinecke's avatar
Martin Reinecke committed
368
        return self.domain.shape
csongor's avatar
csongor committed
369

370
371
    @property
    def dim(self):
Theo Steininger's avatar
Theo Steininger committed
372
        """ Returns the total number of pixel-dimensions the field has.
Theo Steininger's avatar
Theo Steininger committed
373

Theo Steininger's avatar
Theo Steininger committed
374
        Effectively, all values from shape are multiplied.
Theo Steininger's avatar
Theo Steininger committed
375

376
377
378
379
380
        Returns
        -------
        out : int
            The dimension of the Field.
        """
Martin Reinecke's avatar
Martin Reinecke committed
381
        return self.domain.dim
csongor's avatar
csongor committed
382

Theo Steininger's avatar
Theo Steininger committed
383
384
    @property
    def total_volume(self):
Theo Steininger's avatar
Theo Steininger committed
385
386
        """ Returns the total volume of all spaces in the domain.
        """
Martin Reinecke's avatar
Martin Reinecke committed
387
        if len(self.domain) == 0:
Theo Steininger's avatar
Theo Steininger committed
388
            return 0.
Martin Reinecke's avatar
Martin Reinecke committed
389
390
        volume_tuple = tuple(sp.total_volume for sp in self.domain)
        return reduce(lambda x, y: x * y, volume_tuple)
391

Theo Steininger's avatar
Theo Steininger committed
392
393
394
395
    @property
    def real(self):
        """ The real part of the field (data is not copied).
        """
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
396
        return Field(self.domain, self.val.real)
Theo Steininger's avatar
Theo Steininger committed
397
398
399
400
401

    @property
    def imag(self):
        """ The imaginary part of the field (data is not copied).
        """
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
402
        return Field(self.domain, self.val.imag)
Theo Steininger's avatar
Theo Steininger committed
403

Theo Steininger's avatar
Theo Steininger committed
404
    # ---Special unary/binary operations---
405

Martin Reinecke's avatar
Martin Reinecke committed
406
    def copy(self):
407
        """ Returns a full copy of the Field.
Theo Steininger's avatar
Theo Steininger committed
408

Martin Reinecke's avatar
Martin Reinecke committed
409
        The returned object will be an identical copy of the original Field.
Theo Steininger's avatar
Theo Steininger committed
410

411
412
413
414
415
        Returns
        -------
        out : Field
            The output object. An identical copy of 'self'.
        """
Martin Reinecke's avatar
Martin Reinecke committed
416
        return Field(val=self, copy=True)
csongor's avatar
csongor committed
417

418
419
420
421
422
423
    def scalar_weight(self, spaces=None):
        if np.isscalar(spaces):
            return self.domain[spaces].scalar_weight()

        if spaces is None:
            spaces = range(len(self.domain))
Martin Reinecke's avatar
Martin Reinecke committed
424
        res = 1.
425
426
427
428
429
430
431
        for i in spaces:
            tmp = self.domain[i].scalar_weight()
            if tmp is None:
                return None
            res *= tmp
        return res

Theo Steininger's avatar
Theo Steininger committed
432
    def weight(self, power=1, inplace=False, spaces=None):
Theo Steininger's avatar
Theo Steininger committed
433
        """ Weights the pixels of `self` with their invidual pixel-volume.
434
435
436
437

        Parameters
        ----------
        power : number
Theo Steininger's avatar
Theo Steininger committed
438
            The pixels get weighted with the volume-factor**power.
Theo Steininger's avatar
Theo Steininger committed
439

440
        inplace : boolean
Theo Steininger's avatar
Theo Steininger committed
441
442
            If True, `self` will be weighted and returned. Otherwise, a copy
            is made.
Theo Steininger's avatar
Theo Steininger committed
443

Theo Steininger's avatar
Theo Steininger committed
444
445
        spaces : tuple of ints
            Determines on which subspace the operation takes place.
Theo Steininger's avatar
Theo Steininger committed
446

447
448
449
        Returns
        -------
        out : Field
Theo Steininger's avatar
Theo Steininger committed
450
            The weighted field.
451
452

        """
453
        new_field = self if inplace else self.copy()
csongor's avatar
csongor committed
454

csongor's avatar
csongor committed
455
        if spaces is None:
Martin Reinecke's avatar
Martin Reinecke committed
456
            spaces = range(len(self.domain))
Martin Reinecke's avatar
Martin Reinecke committed
457
458
        else:
            spaces = utilities.cast_iseq_to_tuple(spaces)
csongor's avatar
csongor committed
459

460
461
462
463
464
465
        fct = 1.
        for ind in spaces:
            wgt = self.domain[ind].weight()
            if np.isscalar(wgt):
                fct *= wgt
            else:
466
                new_shape = dobj.ones(len(self.shape), dtype=np.int)
Martin Reinecke's avatar
Martin Reinecke committed
467
468
                new_shape[self.domain.axes[ind][0]:
                          self.domain.axes[ind][-1]+1] = wgt.shape
469
                wgt = wgt.reshape(new_shape)
470
                new_field *= wgt**power
471
        fct = fct**power
Martin Reinecke's avatar
Martin Reinecke committed
472
        if fct != 1.:
473
            new_field *= fct
474

475
        return new_field
csongor's avatar
csongor committed
476

Martin Reinecke's avatar
Martin Reinecke committed
477
    def vdot(self, x=None, spaces=None):
Theo Steininger's avatar
Theo Steininger committed
478
        """ Computes the volume-factor-aware dot product of 'self' with x.
Theo Steininger's avatar
Theo Steininger committed
479

480
481
482
        Parameters
        ----------
        x : Field
Theo Steininger's avatar
Theo Steininger committed
483
            The domain of x must contain `self.domain`
Theo Steininger's avatar
Theo Steininger committed
484

Theo Steininger's avatar
Theo Steininger committed
485
486
487
        spaces : tuple of ints
            If the domain of `self` and `x` are not the same, `spaces` specfies
            the mapping.
Theo Steininger's avatar
Theo Steininger committed
488

489
490
491
        Returns
        -------
        out : float, complex
Theo Steininger's avatar
Theo Steininger committed
492

493
        """
494
495
496
        if not isinstance(x, Field):
            raise ValueError("The dot-partner must be an instance of " +
                             "the NIFTy field class")
Theo Steininger's avatar
Theo Steininger committed
497

Martin Reinecke's avatar
Martin Reinecke committed
498
        # Compute the dot respecting the fact of discrete/continuous spaces
499
        fct = 1.
Martin Reinecke's avatar
Martin Reinecke committed
500
501
502
        tmp = self.scalar_weight(spaces)
        if tmp is None:
            y = self.weight(power=1)
503
        else:
Martin Reinecke's avatar
Martin Reinecke committed
504
505
            y = self
            fct = tmp
Theo Steininger's avatar
Theo Steininger committed
506

507
        if spaces is None:
508
            return fct*dobj.vdot(y.val.ravel(), x.val.ravel())
509
510
511
        else:
            # create a diagonal operator which is capable of taking care of the
            # axes-matching
512
            from .operators.diagonal_operator import DiagonalOperator
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
513
514
            diag = DiagonalOperator(y.domain, y.conjugate(), copy=False)
            dotted = diag(x, spaces=spaces)
515
            return fct*dotted.sum(spaces=spaces)
Theo Steininger's avatar
Theo Steininger committed
516

Theo Steininger's avatar
Theo Steininger committed
517
    def norm(self):
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
518
        """ Computes the L2-norm of the field values.
csongor's avatar
csongor committed
519

Theo Steininger's avatar
Theo Steininger committed
520
521
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
522
        norm : float
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
523
            The L2-norm of the field values.
csongor's avatar
csongor committed
524
525

        """
526
        return dobj.sqrt(dobj.abs(self.vdot(x=self)))
csongor's avatar
csongor committed
527

Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
528
    def conjugate(self):
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
529
        """ Returns the complex conjugate of the field.
Theo Steininger's avatar
Theo Steininger committed
530

531
532
533
534
        Returns
        -------
        cc : field
            The complex conjugated field.
csongor's avatar
csongor committed
535
536

        """
Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
537
        return Field(self.domain, self.val.conjugate(), self.dtype)
csongor's avatar
csongor committed
538

Theo Steininger's avatar
Theo Steininger committed
539
    # ---General unary/contraction methods---
540

Theo Steininger's avatar
Theo Steininger committed
541
542
    def __pos__(self):
        return self.copy()
543

Theo Steininger's avatar
Theo Steininger committed
544
    def __neg__(self):
Martin Reinecke's avatar
PEP8    
Martin Reinecke committed
545
        return Field(self.domain, -self.val, self.dtype)
csongor's avatar
csongor committed
546

Theo Steininger's avatar
Theo Steininger committed
547
    def __abs__(self):
548
        return Field(self.domain, dobj.abs(self.val), self.dtype)
csongor's avatar
csongor committed
549

550
    def _contraction_helper(self, op, spaces):
Theo Steininger's avatar
Theo Steininger committed
551
        if spaces is None:
552
            return getattr(self.val, op)()
Martin Reinecke's avatar
Martin Reinecke committed
553
554
        else:
            spaces = utilities.cast_iseq_to_tuple(spaces)
csongor's avatar
csongor committed
555

Martin Reinecke's avatar
Martin Reinecke committed
556
        axes_list = tuple(self.domain.axes[sp_index] for sp_index in spaces)
557

Martin Reinecke's avatar
Martin Reinecke committed
558
        if len(axes_list) > 0:
Theo Steininger's avatar
Theo Steininger committed
559
            axes_list = reduce(lambda x, y: x+y, axes_list)
csongor's avatar
csongor committed
560

Martin Reinecke's avatar
stage1    
Martin Reinecke committed
561
        # perform the contraction on the data
562
        data = getattr(self.val, op)(axis=axes_list)
csongor's avatar
csongor committed
563

Theo Steininger's avatar
Theo Steininger committed
564
565
566
        # check if the result is scalar or if a result_field must be constr.
        if np.isscalar(data):
            return data
csongor's avatar
csongor committed
567
        else:
Theo Steininger's avatar
Theo Steininger committed
568
            return_domain = tuple(self.domain[i]
Martin Reinecke's avatar
Martin Reinecke committed
569
                                  for i in range(len(self.domain))
Theo Steininger's avatar
Theo Steininger committed
570
                                  if i not in spaces)
571

Martin Reinecke's avatar
updates    
Martin Reinecke committed
572
            return Field(domain=return_domain, val=data, copy=False)
csongor's avatar
csongor committed
573

574
575
    def sum(self, spaces=None):
        return self._contraction_helper('sum', spaces)
csongor's avatar
csongor committed
576

577
578
    def prod(self, spaces=None):
        return self._contraction_helper('prod', spaces)
csongor's avatar
csongor committed
579

580
581
    def all(self, spaces=None):
        return self._contraction_helper('all', spaces)
csongor's avatar
csongor committed
582

583
584
    def any(self, spaces=None):
        return self._contraction_helper('any', spaces)
csongor's avatar
csongor committed
585

586
587
    def min(self, spaces=None):
        return self._contraction_helper('min', spaces)
csongor's avatar
csongor committed
588

589
590
    def max(self, spaces=None):
        return self._contraction_helper('max', spaces)
csongor's avatar
csongor committed
591

592
593
    def mean(self, spaces=None):
        return self._contraction_helper('mean', spaces)
csongor's avatar
csongor committed
594

595
596
    def var(self, spaces=None):
        return self._contraction_helper('var', spaces)
csongor's avatar
csongor committed
597

598
599
    def std(self, spaces=None):
        return self._contraction_helper('std', spaces)
csongor's avatar
csongor committed
600

601
602
603
604
605
606
607
    def copy_content_from(self, other):
        if not isinstance(other, Field):
            raise TypeError("argument must be a Field")
        if other.domain != self.domain:
            raise ValueError("domains are incompatible.")
        self.val[()] = other.val

Theo Steininger's avatar
Theo Steininger committed
608
    # ---General binary methods---
csongor's avatar
csongor committed
609

610
    def _binary_helper(self, other, op):
csongor's avatar
csongor committed
611
        # if other is a field, make sure that the domains match
612
        if isinstance(other, Field):
613
614
            if other.domain != self.domain:
                raise ValueError("domains are incompatible.")
Martin Reinecke's avatar
Martin Reinecke committed
615
616
            tval = getattr(self.val, op)(other.val)
            return self if tval is self.val else Field(self.domain, tval)
csongor's avatar
csongor committed
617

Martin Reinecke's avatar
Martin Reinecke committed
618
619
        tval = getattr(self.val, op)(other)
        return self if tval is self.val else Field(self.domain, tval)
csongor's avatar
csongor committed
620
621

    def __add__(self, other):
Theo Steininger's avatar
Theo Steininger committed
622
        return self._binary_helper(other, op='__add__')
623

624
    def __radd__(self, other):
Theo Steininger's avatar
Theo Steininger committed
625
        return self._binary_helper(other, op='__radd__')
csongor's avatar
csongor committed
626
627

    def __iadd__(self, other):
628
        return self._binary_helper(other, op='__iadd__')
csongor's avatar
csongor committed
629
630

    def __sub__(self, other):
Theo Steininger's avatar
Theo Steininger committed
631
        return self._binary_helper(other, op='__sub__')
csongor's avatar
csongor committed
632
633

    def __rsub__(self, other):
Theo Steininger's avatar
Theo Steininger committed
634
        return self._binary_helper(other, op='__rsub__')
csongor's avatar
csongor committed
635
636

    def __isub__(self, other):
637
        return self._binary_helper(other, op='__isub__')
csongor's avatar
csongor committed
638
639

    def __mul__(self, other):
Theo Steininger's avatar
Theo Steininger committed
640
        return self._binary_helper(other, op='__mul__')
641

642
    def __rmul__(self, other):
Theo Steininger's avatar
Theo Steininger committed
643
        return self._binary_helper(other, op='__rmul__')
csongor's avatar
csongor committed
644
645

    def __imul__(self, other):
646
        return self._binary_helper(other, op='__imul__')
csongor's avatar
csongor committed
647
648

    def __div__(self, other):
Theo Steininger's avatar
Theo Steininger committed
649
        return self._binary_helper(other, op='__div__')
csongor's avatar
csongor committed
650

Martin Reinecke's avatar
Martin Reinecke committed
651
652
653
    def __truediv__(self, other):
        return self._binary_helper(other, op='__truediv__')

csongor's avatar
csongor committed
654
    def __rdiv__(self, other):
Theo Steininger's avatar
Theo Steininger committed
655
        return self._binary_helper(other, op='__rdiv__')
csongor's avatar
csongor committed
656

Martin Reinecke's avatar
Martin Reinecke committed
657
658
659
    def __rtruediv__(self, other):
        return self._binary_helper(other, op='__rtruediv__')

csongor's avatar
csongor committed
660
    def __idiv__(self, other):
661
        return self._binary_helper(other, op='__idiv__')
662

csongor's avatar
csongor committed
663
    def __pow__(self, other):
Theo Steininger's avatar
Theo Steininger committed
664
        return self._binary_helper(other, op='__pow__')
csongor's avatar
csongor committed
665
666

    def __rpow__(self, other):
Theo Steininger's avatar
Theo Steininger committed
667
        return self._binary_helper(other, op='__rpow__')
csongor's avatar
csongor committed
668
669

    def __ipow__(self, other):
670
        return self._binary_helper(other, op='__ipow__')
csongor's avatar
csongor committed
671
672

    def __lt__(self, other):
Theo Steininger's avatar
Theo Steininger committed
673
        return self._binary_helper(other, op='__lt__')
csongor's avatar
csongor committed
674
675

    def __le__(self, other):
Theo Steininger's avatar
Theo Steininger committed
676
        return self._binary_helper(other, op='__le__')
csongor's avatar
csongor committed
677
678
679
680
681

    def __ne__(self, other):
        if other is None:
            return True
        else:
Theo Steininger's avatar
Theo Steininger committed
682
            return self._binary_helper(other, op='__ne__')
csongor's avatar
csongor committed
683
684
685
686
687

    def __eq__(self, other):
        if other is None:
            return False
        else:
Theo Steininger's avatar
Theo Steininger committed
688
            return self._binary_helper(other, op='__eq__')
csongor's avatar
csongor committed
689
690

    def __ge__(self, other):
Theo Steininger's avatar
Theo Steininger committed
691
        return self._binary_helper(other, op='__ge__')
csongor's avatar
csongor committed
692
693

    def __gt__(self, other):
Theo Steininger's avatar
Theo Steininger committed
694
695
696
        return self._binary_helper(other, op='__gt__')

    def __repr__(self):
Martin Reinecke's avatar
Martin Reinecke committed
697
        return "<nifty2go.Field>"
Theo Steininger's avatar
Theo Steininger committed
698
699
700
701

    def __str__(self):
        minmax = [self.min(), self.max()]
        mean = self.mean()
Martin Reinecke's avatar
Martin Reinecke committed
702
        return "nifty2go.Field instance\n- domain      = " + \
Theo Steininger's avatar
Theo Steininger committed
703
               repr(self.domain) + \
704
               "\n- val         = " + repr(self.val) + \
Theo Steininger's avatar
Theo Steininger committed
705
706
               "\n  - min.,max. = " + str(minmax) + \
               "\n  - mean = " + str(mean)