field.py 21.5 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
#
Martin Reinecke's avatar
Martin Reinecke committed
14
# Copyright(C) 2013-2018 Max-Planck-Society
Theo Steininger's avatar
Theo Steininger committed
15
16
17
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
# and financially supported by the Studienstiftung des deutschen Volkes.
18

19
from __future__ import absolute_import, division, print_function
20

csongor's avatar
csongor committed
21
import numpy as np
22
23
24

from . import dobj, utilities
from .compat import *
Martin Reinecke's avatar
Martin Reinecke committed
25
from .domain_tuple import DomainTuple
26

27

Martin Reinecke's avatar
Martin Reinecke committed
28
class Field(object):
Martin Reinecke's avatar
Martin Reinecke committed
29
30
    _scalar_dom = DomainTuple.scalar_domain()

Theo Steininger's avatar
Theo Steininger committed
31
32
    """ The discrete representation of a continuous field over multiple spaces.

Martin Reinecke's avatar
Martin Reinecke committed
33
    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.
Theo Steininger's avatar
Theo Steininger committed
35

36
37
    Parameters
    ----------
38
39
    domain : DomainTuple
        the domain of the new Field
Theo Steininger's avatar
Theo Steininger committed
40

41
42
43
    val : data_object
        This object's global shape must match the domain shape
        After construction, the object will no longer be writeable!
Martin Reinecke's avatar
Martin Reinecke committed
44
45
46
47

    Notes
    -----
    If possible, do not invoke the constructor directly, but use one of the
48
    many convenience functions for Field construction!
49
    """
50

51
52
53
    def __init__(self, domain, val):
        if not isinstance(domain, DomainTuple):
            raise TypeError("domain must be of type DomainTuple")
Martin Reinecke's avatar
Martin Reinecke committed
54
        if type(val) is not dobj.data_object:
Martin Reinecke's avatar
Martin Reinecke committed
55
            if np.isscalar(val):
Martin Reinecke's avatar
Martin Reinecke committed
56
                val = dobj.full(domain.shape, val)
Martin Reinecke's avatar
Martin Reinecke committed
57
58
            else:
                raise TypeError("val must be of type dobj.data_object")
59
        if domain.shape != val.shape:
Martin Reinecke's avatar
Martin Reinecke committed
60
            raise ValueError("shape mismatch between val and domain")
61
62
        self._domain = domain
        self._val = val
63
        dobj.lock(self._val)
Martin Reinecke's avatar
Martin Reinecke committed
64

Martin Reinecke's avatar
Martin Reinecke committed
65
66
67
68
    @staticmethod
    def scalar(val):
        return Field(Field._scalar_dom, val)

69
70
71
72
73
74
75
    # prevent implicit conversion to bool
    def __nonzero__(self):
        raise TypeError("Field does not support implicit conversion to bool")

    def __bool__(self):
        raise TypeError("Field does not support implicit conversion to bool")

76
    @staticmethod
77
    def full(domain, val):
Martin Reinecke's avatar
Martin Reinecke committed
78
79
80
81
82
83
84
85
86
87
88
89
90
91
        """Creates a Field with a given domain, filled with a constant value.

        Parameters
        ----------
        domain : Domain, tuple of Domain, or DomainTuple
            domain of the new Field
        val : float/complex/int scalar
            fill value. Data type of the field is inferred from val.

        Returns
        -------
        Field
            the newly created field
        """
92
93
        if not np.isscalar(val):
            raise TypeError("val must be a scalar")
94
95
96
        if not (np.isreal(val) or np.iscomplex(val)):
            raise TypeError("need arithmetic scalar")
        domain = DomainTuple.make(domain)
Martin Reinecke's avatar
Martin Reinecke committed
97
        return Field(domain, val)
98

99
    @staticmethod
100
    def from_global_data(domain, arr, sum_up=False):
Martin Reinecke's avatar
Martin Reinecke committed
101
102
103
104
105
106
107
108
109
        """Returns a Field constructed from `domain` and `arr`.

        Parameters
        ----------
        domain : DomainTuple, tuple of Domain, or Domain
            the domain of the new Field
        arr : numpy.ndarray
            The data content to be used for the new Field.
            Its shape must match the shape of `domain`.
110
111
112
113
114
        sum_up : bool, optional
            If True, the contents of `arr` are summed up over all MPI tasks
            (if any), and the sum is used as data content.
            If False, the contens of `arr` are used directly, and must be
            identical on all MPI tasks.
Martin Reinecke's avatar
Martin Reinecke committed
115
        """
116
117
        return Field(DomainTuple.make(domain),
                     dobj.from_global_data(arr, sum_up))
118

Martin Reinecke's avatar
Martin Reinecke committed
119
120
    @staticmethod
    def from_local_data(domain, arr):
121
        return Field(DomainTuple.make(domain),
Martin Reinecke's avatar
Martin Reinecke committed
122
                     dobj.from_local_data(domain.shape, arr))
Martin Reinecke's avatar
Martin Reinecke committed
123

124
    def to_global_data(self):
Martin Reinecke's avatar
Martin Reinecke committed
125
126
127
128
129
130
131
        """Returns an array containing the full data of the field.

        Returns
        -------
        numpy.ndarray : array containing all field entries.
            Its shape is identical to `self.shape`.
        """
132
133
        return dobj.to_global_data(self._val)

134
135
136
137
138
139
140
141
142
143
    def to_global_data_rw(self):
        """Returns a modifiable array containing the full data of the field.

        Returns
        -------
        numpy.ndarray : array containing all field entries, which can be
            modified. Its shape is identical to `self.shape`.
        """
        return dobj.to_global_data_rw(self._val)

Martin Reinecke's avatar
Martin Reinecke committed
144
145
    @property
    def local_data(self):
Martin Reinecke's avatar
Martin Reinecke committed
146
147
148
149
150
        """numpy.ndarray : locally residing field data

        Returns a handle to the part of the array data residing on the local
        task (or to the entore array if MPI is not active).
        """
Martin Reinecke's avatar
Martin Reinecke committed
151
152
        return dobj.local_data(self._val)

Martin Reinecke's avatar
Martin Reinecke committed
153
    def cast_domain(self, new_domain):
Martin Reinecke's avatar
Martin Reinecke committed
154
155
156
157
158
159
160
161
162
163
164
165
166
        """Returns a field with the same data, but a different domain

        Parameters
        ----------
        new_domain : Domain, tuple of Domain, or DomainTuple
            The domain for the returned field. Must be shape-compatible to
            `self`.

        Returns
        -------
        Field
            Field living on `new_domain`, but with the same data as `self`.
        """
167
        return Field(DomainTuple.make(new_domain), self._val)
168

Martin Reinecke's avatar
Martin Reinecke committed
169
170
    @staticmethod
    def from_random(random_type, domain, dtype=np.float64, **kwargs):
171
172
173
174
        """ Draws a random field with the given parameters.

        Parameters
        ----------
Martin Reinecke's avatar
Martin Reinecke committed
175
176
        random_type : 'pm1', 'normal', or 'uniform'
            The random distribution to use.
Martin Reinecke's avatar
Martin Reinecke committed
177
        domain : DomainTuple
178
179
180
            The domain of the output random field
        dtype : type
            The datatype of the output random field
Theo Steininger's avatar
Theo Steininger committed
181

182
183
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
184
        Field
Martin Reinecke's avatar
Martin Reinecke committed
185
            The newly created Field.
186
        """
Martin Reinecke's avatar
Martin Reinecke committed
187
        domain = DomainTuple.make(domain)
Martin Reinecke's avatar
Martin Reinecke committed
188
189
190
        return Field(domain=domain,
                     val=dobj.from_random(random_type, dtype=dtype,
                                          shape=domain.shape, **kwargs))
191

Theo Steininger's avatar
Theo Steininger committed
192
193
    @property
    def val(self):
Martin Reinecke's avatar
Martin Reinecke committed
194
195
196
197
198
199
        """dobj.data_object : the data object storing the field's entries

        Notes
        -----
        This property is intended for low-level, internal use only. Do not use
        from outside of NIFTy's core; there should be better alternatives.
200
        """
Martin Reinecke's avatar
Martin Reinecke committed
201
        return self._val
csongor's avatar
csongor committed
202

Martin Reinecke's avatar
Martin Reinecke committed
203
204
    @property
    def dtype(self):
Martin Reinecke's avatar
Martin Reinecke committed
205
        """type : the data type of the field's entries"""
Martin Reinecke's avatar
Martin Reinecke committed
206
207
        return self._val.dtype

Martin Reinecke's avatar
Martin Reinecke committed
208
209
    @property
    def domain(self):
Martin Reinecke's avatar
Martin Reinecke committed
210
        """DomainTuple : the field's domain"""
Martin Reinecke's avatar
Martin Reinecke committed
211
212
        return self._domain

213
214
    @property
    def shape(self):
Martin Reinecke's avatar
Martin Reinecke committed
215
        """tuple of int : the concatenated shapes of all sub-domains"""
Martin Reinecke's avatar
Martin Reinecke committed
216
        return self._domain.shape
csongor's avatar
csongor committed
217

218
    @property
Martin Reinecke's avatar
Martin Reinecke committed
219
    def size(self):
Martin Reinecke's avatar
Martin Reinecke committed
220
        """int : total number of pixels in the field"""
Martin Reinecke's avatar
Martin Reinecke committed
221
        return self._domain.size
csongor's avatar
csongor committed
222

Theo Steininger's avatar
Theo Steininger committed
223
224
    @property
    def real(self):
Martin Reinecke's avatar
Martin Reinecke committed
225
        """Field : The real part of the field"""
Martin Reinecke's avatar
Martin Reinecke committed
226
227
228
        if utilities.iscomplextype(self.dtype):
            return Field(self._domain, self._val.real)
        return self
Theo Steininger's avatar
Theo Steininger committed
229
230
231

    @property
    def imag(self):
Martin Reinecke's avatar
Martin Reinecke committed
232
        """Field : The imaginary part of the field"""
Martin Reinecke's avatar
Martin Reinecke committed
233
        if not utilities.iscomplextype(self.dtype):
234
            raise ValueError(".imag called on a non-complex Field")
Martin Reinecke's avatar
Martin Reinecke committed
235
        return Field(self._domain, self._val.imag)
Theo Steininger's avatar
Theo Steininger committed
236

237
    def scalar_weight(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
238
239
240
241
242
243
244
245
246
247
248
249
250
251
        """Returns the uniform volume element for a sub-domain of `self`.

        Parameters
        ----------
        spaces : int, tuple of int or None
            indices of the sub-domains of the field's domain to be considered.
            If `None`, the entire domain is used.

        Returns
        -------
        float or None
            if the requested sub-domain has a uniform volume element, it is
            returned. Otherwise, `None` is returned.
        """
252
        if np.isscalar(spaces):
Martin Reinecke's avatar
Martin Reinecke committed
253
            return self._domain[spaces].scalar_dvol
254
255

        if spaces is None:
Martin Reinecke's avatar
Martin Reinecke committed
256
            spaces = range(len(self._domain))
Martin Reinecke's avatar
Martin Reinecke committed
257
        res = 1.
258
        for i in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
259
            tmp = self._domain[i].scalar_dvol
260
261
262
263
264
            if tmp is None:
                return None
            res *= tmp
        return res

Martin Reinecke's avatar
Martin Reinecke committed
265
    def total_volume(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
266
267
268
269
270
271
272
273
274
275
276
277
278
        """Returns the total volume of a sub-domain of `self`.

        Parameters
        ----------
        spaces : int, tuple of int or None
            indices of the sub-domains of the field's domain to be considered.
            If `None`, the entire domain is used.

        Returns
        -------
        float
            the total volume of the requested sub-domain.
        """
Martin Reinecke's avatar
Martin Reinecke committed
279
        if np.isscalar(spaces):
Martin Reinecke's avatar
Martin Reinecke committed
280
            return self._domain[spaces].total_volume
Martin Reinecke's avatar
Martin Reinecke committed
281
282
283
284
285

        if spaces is None:
            spaces = range(len(self._domain))
        res = 1.
        for i in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
286
            res *= self._domain[i].total_volume
Martin Reinecke's avatar
Martin Reinecke committed
287
288
        return res

289
    def weight(self, power=1, spaces=None):
Theo Steininger's avatar
Theo Steininger committed
290
        """ Weights the pixels of `self` with their invidual pixel-volume.
291
292
293
294

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

Martin Reinecke's avatar
Martin Reinecke committed
297
298
299
        spaces : None, int or tuple of int
            Determines on which sub-domain the operation takes place.
            If None, the entire domain is used.
Theo Steininger's avatar
Theo Steininger committed
300

301
302
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
303
        Field
Theo Steininger's avatar
Theo Steininger committed
304
            The weighted field.
305
        """
306
        aout = self.local_data.copy()
csongor's avatar
csongor committed
307

Martin Reinecke's avatar
Martin Reinecke committed
308
        spaces = utilities.parse_spaces(spaces, len(self._domain))
csongor's avatar
csongor committed
309

310
311
        fct = 1.
        for ind in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
312
            wgt = self._domain[ind].dvol
313
314
315
            if np.isscalar(wgt):
                fct *= wgt
            else:
Martin Reinecke's avatar
Martin Reinecke committed
316
                new_shape = np.ones(len(self.shape), dtype=np.int)
Martin Reinecke's avatar
Martin Reinecke committed
317
318
                new_shape[self._domain.axes[ind][0]:
                          self._domain.axes[ind][-1]+1] = wgt.shape
319
                wgt = wgt.reshape(new_shape)
Martin Reinecke's avatar
Martin Reinecke committed
320
321
                if dobj.distaxis(self._val) >= 0 and ind == 0:
                    # we need to distribute the weights along axis 0
Martin Reinecke's avatar
fixes  
Martin Reinecke committed
322
                    wgt = dobj.local_data(dobj.from_global_data(wgt))
323
                aout *= wgt**power
324
        fct = fct**power
Martin Reinecke's avatar
Martin Reinecke committed
325
        if fct != 1.:
326
            aout *= fct
327

328
        return Field.from_local_data(self._domain, aout)
csongor's avatar
csongor committed
329

Martin Reinecke's avatar
Martin Reinecke committed
330
    def vdot(self, x=None, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
331
        """ Computes the dot product of 'self' with x.
Theo Steininger's avatar
Theo Steininger committed
332

333
334
335
        Parameters
        ----------
        x : Field
336
            x must live on the same domain as `self`.
Theo Steininger's avatar
Theo Steininger committed
337

Martin Reinecke's avatar
Martin Reinecke committed
338
        spaces : None, int or tuple of int (default: None)
339
340
            The dot product is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.
Theo Steininger's avatar
Theo Steininger committed
341

342
343
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
344
        float, complex, either scalar (for full dot products)
345
                              or Field (for partial dot products)
346
        """
347
        if not isinstance(x, Field):
348
349
            raise TypeError("The dot-partner must be an instance of " +
                            "the NIFTy field class")
Theo Steininger's avatar
Theo Steininger committed
350

351
        if x._domain is not self._domain:
352
            raise ValueError("Domain mismatch")
Theo Steininger's avatar
Theo Steininger committed
353

Martin Reinecke's avatar
Martin Reinecke committed
354
        ndom = len(self._domain)
355
356
357
        spaces = utilities.parse_spaces(spaces, ndom)

        if len(spaces) == ndom:
Martin Reinecke's avatar
Martin Reinecke committed
358
            return dobj.vdot(self._val, x._val)
359
360
        # If we arrive here, we have to do a partial dot product.
        # For the moment, do this the explicit, non-optimized way
Martin Reinecke's avatar
Martin Reinecke committed
361
        return (self.conjugate()*x).sum(spaces=spaces)
Theo Steininger's avatar
Theo Steininger committed
362

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

Theo Steininger's avatar
Theo Steininger committed
366
367
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
368
        float
Martin Reinecke's avatar
tweaks  
Martin Reinecke committed
369
            The L2-norm of the field values.
csongor's avatar
csongor committed
370
        """
371
        return np.sqrt(abs(self.vdot(x=self)))
csongor's avatar
csongor committed
372

373
374
375
376
377
378
379
380
381
382
    def squared_norm(self):
        """ Computes the square of the L2-norm of the field values.

        Returns
        -------
        float
            The square of the L2-norm of the field values.
        """
        return abs(self.vdot(x=self))

Martin Reinecke's avatar
tweaks  
Martin Reinecke committed
383
    def conjugate(self):
Martin Reinecke's avatar
cleanup  
Martin Reinecke committed
384
        """ Returns the complex conjugate of the field.
Theo Steininger's avatar
Theo Steininger committed
385

386
387
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
388
389
        Field
            The complex conjugated field.
csongor's avatar
csongor committed
390
        """
Martin Reinecke's avatar
Martin Reinecke committed
391
        if utilities.iscomplextype(self._val.dtype):
Martin Reinecke's avatar
Martin Reinecke committed
392
393
            return Field(self._domain, self._val.conjugate())
        return self
csongor's avatar
csongor committed
394

Theo Steininger's avatar
Theo Steininger committed
395
    # ---General unary/contraction methods---
396

Theo Steininger's avatar
Theo Steininger committed
397
    def __pos__(self):
398
        return self
399

Theo Steininger's avatar
Theo Steininger committed
400
    def __neg__(self):
Martin Reinecke's avatar
Martin Reinecke committed
401
        return Field(self._domain, -self._val)
csongor's avatar
csongor committed
402

Theo Steininger's avatar
Theo Steininger committed
403
    def __abs__(self):
Martin Reinecke's avatar
Martin Reinecke committed
404
        return Field(self._domain, abs(self._val))
csongor's avatar
csongor committed
405

406
    def _contraction_helper(self, op, spaces):
Theo Steininger's avatar
Theo Steininger committed
407
        if spaces is None:
Martin Reinecke's avatar
Martin Reinecke committed
408
            return getattr(self._val, op)()
409

Martin Reinecke's avatar
Martin Reinecke committed
410
        spaces = utilities.parse_spaces(spaces, len(self._domain))
csongor's avatar
csongor committed
411

Martin Reinecke's avatar
Martin Reinecke committed
412
        axes_list = tuple(self._domain.axes[sp_index] for sp_index in spaces)
413

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

Martin Reinecke's avatar
stage1  
Martin Reinecke committed
417
        # perform the contraction on the data
Martin Reinecke's avatar
Martin Reinecke committed
418
        data = getattr(self._val, op)(axis=axes_list)
csongor's avatar
csongor committed
419

Theo Steininger's avatar
Theo Steininger committed
420
421
422
        # 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
423
        else:
Martin Reinecke's avatar
Martin Reinecke committed
424
            return_domain = tuple(dom
Martin Reinecke's avatar
Martin Reinecke committed
425
                                  for i, dom in enumerate(self._domain)
Theo Steininger's avatar
Theo Steininger committed
426
                                  if i not in spaces)
427

428
            return Field(DomainTuple.make(return_domain), data)
csongor's avatar
csongor committed
429

430
    def sum(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
431
432
433
434
435
436
437
438
439
440
441
442
443
444
        """Sums up over the sub-domains given by `spaces`.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The summation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the summation. If it is carried out over the entire
            domain, this is a scalar, otherwise a Field.
        """
445
        return self._contraction_helper('sum', spaces)
csongor's avatar
csongor committed
446

447
    def integrate(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
        """Integrates over the sub-domains given by `spaces`.

        Integration is performed by summing over `self` multiplied by its
        volume factors.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The summation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the integration. If it is carried out over the
            entire domain, this is a scalar, otherwise a Field.
        """
Martin Reinecke's avatar
Martin Reinecke committed
465
466
467
468
469
        swgt = self.scalar_weight(spaces)
        if swgt is not None:
            res = self.sum(spaces)
            res *= swgt
            return res
470
471
472
        tmp = self.weight(1, spaces=spaces)
        return tmp.sum(spaces)

473
    def prod(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
474
475
476
477
478
479
480
481
482
483
484
485
486
487
        """Computes the product over the sub-domains given by `spaces`.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the product. If it is carried out over the entire
            domain, this is a scalar, otherwise a Field.
        """
488
        return self._contraction_helper('prod', spaces)
csongor's avatar
csongor committed
489

490
491
    def all(self, spaces=None):
        return self._contraction_helper('all', spaces)
csongor's avatar
csongor committed
492

493
494
    def any(self, spaces=None):
        return self._contraction_helper('any', spaces)
csongor's avatar
csongor committed
495

496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
#     def min(self, spaces=None):
#         """Determines the minimum over the sub-domains given by `spaces`.
#
#         Parameters
#         ----------
#         spaces : None, int or tuple of int (default: None)
#             The operation is only carried out over the sub-domains in this
#             tuple. If None, it is carried out over all sub-domains.
#
#         Returns
#         -------
#         Field or scalar
#             The result of the operation. If it is carried out over the entire
#             domain, this is a scalar, otherwise a Field.
#         """
#         return self._contraction_helper('min', spaces)
#
#     def max(self, spaces=None):
#         """Determines the maximum over the sub-domains given by `spaces`.
#
#         Parameters
#         ----------
#         spaces : None, int or tuple of int (default: None)
#             The operation is only carried out over the sub-domains in this
#             tuple. If None, it is carried out over all sub-domains.
#
#         Returns
#         -------
#         Field or scalar
#             The result of the operation. If it is carried out over the entire
#             domain, this is a scalar, otherwise a Field.
#         """
#         return self._contraction_helper('max', spaces)
csongor's avatar
csongor committed
529

530
    def mean(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
        """Determines the mean over the sub-domains given by `spaces`.

        ``x.mean(spaces)`` is equivalent to
        ``x.integrate(spaces)/x.total_volume(spaces)``.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the operation. If it is carried out over the entire
            domain, this is a scalar, otherwise a Field.
        """
548
549
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('mean', spaces)
Martin Reinecke's avatar
Martin Reinecke committed
550
        # MR FIXME: not very efficient
551
552
        # MR FIXME: do we need "spaces" here?
        tmp = self.weight(1, spaces)
Martin Reinecke's avatar
Martin Reinecke committed
553
        return tmp.sum(spaces)*(1./tmp.total_volume(spaces))
csongor's avatar
csongor committed
554

555
    def var(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
556
557
558
559
560
561
562
563
564
565
566
567
568
569
        """Determines the variance over the sub-domains given by `spaces`.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the operation. If it is carried out over the entire
            domain, this is a scalar, otherwise a Field.
        """
570
571
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('var', spaces)
Martin Reinecke's avatar
Martin Reinecke committed
572
573
        # MR FIXME: not very efficient or accurate
        m1 = self.mean(spaces)
Martin Reinecke's avatar
Martin Reinecke committed
574
        if utilities.iscomplextype(self.dtype):
575
            sq = abs(self-m1)**2
Martin Reinecke's avatar
Martin Reinecke committed
576
        else:
577
578
            sq = (self-m1)**2
        return sq.mean(spaces)
csongor's avatar
csongor committed
579

580
    def std(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
        """Determines the standard deviation over the sub-domains given by
        `spaces`.

        ``x.std(spaces)`` is equivalent to ``sqrt(x.var(spaces))``.

        Parameters
        ----------
        spaces : None, int or tuple of int (default: None)
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
        Field or scalar
            The result of the operation. If it is carried out over the entire
            domain, this is a scalar, otherwise a Field.
        """
598
        from .sugar import sqrt
599
600
601
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('std', spaces)
        return sqrt(self.var(spaces))
csongor's avatar
csongor committed
602

Theo Steininger's avatar
Theo Steininger committed
603
    def __repr__(self):
Philipp Arras's avatar
Philipp Arras committed
604
        return "<nifty5.Field>"
Theo Steininger's avatar
Theo Steininger committed
605
606

    def __str__(self):
Philipp Arras's avatar
Philipp Arras committed
607
        return "nifty5.Field instance\n- domain      = " + \
608
               self._domain.__str__() + \
Martin Reinecke's avatar
Martin Reinecke committed
609
               "\n- val         = " + repr(self._val)
Martin Reinecke's avatar
cleanup  
Martin Reinecke committed
610

Martin Reinecke's avatar
more  
Martin Reinecke committed
611
612
613
614
615
616
617
    def extract(self, dom):
        if dom is not self._domain:
            raise ValueError("domain mismatch")
        return self

    def unite(self, other):
        return self + other
618

Martin Reinecke's avatar
Martin Reinecke committed
619
    def positive_tanh(self):
Martin Reinecke's avatar
Martin Reinecke committed
620
621
622
623
624
625
        return 0.5*(1.+self.tanh())


for op in ["__add__", "__radd__",
           "__sub__", "__rsub__",
           "__mul__", "__rmul__",
626
627
628
629
           "__div__", "__rdiv__",
           "__truediv__", "__rtruediv__",
           "__floordiv__", "__rfloordiv__",
           "__pow__", "__rpow__",
630
631
632
           "__lt__", "__le__", "__gt__", "__ge__", "__eq__", "__ne__"]:
    def func(op):
        def func2(self, other):
Martin Reinecke's avatar
cleanup  
Martin Reinecke committed
633
634
            # if other is a field, make sure that the domains match
            if isinstance(other, Field):
635
                if other._domain is not self._domain:
Martin Reinecke's avatar
cleanup  
Martin Reinecke committed
636
                    raise ValueError("domains are incompatible.")
Martin Reinecke's avatar
Martin Reinecke committed
637
                tval = getattr(self._val, op)(other._val)
638
                return Field(self._domain, tval)
639
            if np.isscalar(other):
Martin Reinecke's avatar
Martin Reinecke committed
640
                tval = getattr(self._val, op)(other)
641
                return Field(self._domain, tval)
Martin Reinecke's avatar
cleanup  
Martin Reinecke committed
642
            return NotImplemented
643
644
        return func2
    setattr(Field, op, func(op))
645
646
647
648
649
650
651
652
653

for op in ["__iadd__", "__isub__", "__imul__", "__idiv__",
           "__itruediv__", "__ifloordiv__", "__ipow__"]:
    def func(op):
        def func2(self, other):
            raise TypeError(
                "In-place operations are deliberately not supported")
        return func2
    setattr(Field, op, func(op))
Martin Reinecke's avatar
Martin Reinecke committed
654
655
656
657

for f in ["sqrt", "exp", "log", "tanh"]:
    def func(f):
        def func2(self):
Martin Reinecke's avatar
Martin Reinecke committed
658
659
            fu = getattr(dobj, f)
            return Field(domain=self._domain, val=fu(self.val))
Martin Reinecke's avatar
Martin Reinecke committed
660
661
        return func2
    setattr(Field, f, func(f))