1
# Incompatibility in that getmask and a.mask returns nomask
4
from numpy.core.ma import *
5
import numpy.core.ma as nca
1
"""MA: a facility for dealing with missing observations
2
MA is generally used as a numpy.array look-alike.
5
Copyright 1999, 2000, 2001 Regents of the University of California.
6
Released for unlimited redistribution.
7
Adapted for numpy_core 2005 by Travis Oliphant and
13
import numpy.core.umath as umath
14
import numpy.core.fromnumeric as fromnumeric
15
from numpy.core.numeric import newaxis, ndarray, inf
16
from numpy.core.fromnumeric import amax, amin
17
from numpy.core.numerictypes import bool_, typecodes
18
import numpy.core.numeric as numeric
21
# Ufunc domain lookup for __array_wrap__
23
# Ufunc fills lookup for __array__
28
divide_tolerance = 1.e-35
30
class MAError (Exception):
31
def __init__ (self, args=None):
34
# The .args attribute must be a tuple.
35
if not isinstance(args, tuple):
39
"Calculate the string representation"
40
return str(self.args[0])
43
class _MaskedPrintOption:
44
"One instance of this class, masked_print_option, is created."
45
def __init__ (self, display):
46
"Create the masked print option object."
47
self.set_display(display)
51
"Show what prints for masked values."
54
def set_display (self, s):
55
"set_display(s) sets what prints for masked values."
59
"Is the use of the display value enabled?"
62
def enable(self, flag=1):
63
"Set the enabling flag to flag."
67
return str(self._display)
71
#if you single index into a masked location you get this object.
72
masked_print_option = _MaskedPrintOption('--')
74
# Use single element arrays or scalars.
75
default_real_fill_value = 1.e20
76
default_complex_fill_value = 1.e20 + 0.0j
77
default_character_fill_value = '-'
78
default_integer_fill_value = 999999
79
default_object_fill_value = '?'
81
def default_fill_value (obj):
82
"Function to calculate default fill value for an object."
83
if isinstance(obj, types.FloatType):
84
return default_real_fill_value
85
elif isinstance(obj, types.IntType) or isinstance(obj, types.LongType):
86
return default_integer_fill_value
87
elif isinstance(obj, types.StringType):
88
return default_character_fill_value
89
elif isinstance(obj, types.ComplexType):
90
return default_complex_fill_value
91
elif isinstance(obj, MaskedArray) or isinstance(obj, ndarray):
93
if x in typecodes['Float']:
94
return default_real_fill_value
95
if x in typecodes['Integer']:
96
return default_integer_fill_value
97
if x in typecodes['Complex']:
98
return default_complex_fill_value
99
if x in typecodes['Character']:
100
return default_character_fill_value
101
if x in typecodes['UnsignedInteger']:
102
return umath.absolute(default_integer_fill_value)
103
return default_object_fill_value
105
return default_object_fill_value
107
def minimum_fill_value (obj):
108
"Function to calculate default fill value suitable for taking minima."
109
if isinstance(obj, types.FloatType):
111
elif isinstance(obj, types.IntType) or isinstance(obj, types.LongType):
113
elif isinstance(obj, MaskedArray) or isinstance(obj, ndarray):
115
if x in typecodes['Float']:
117
if x in typecodes['Integer']:
119
if x in typecodes['UnsignedInteger']:
122
raise TypeError, 'Unsuitable type for calculating minimum.'
124
def maximum_fill_value (obj):
125
"Function to calculate default fill value suitable for taking maxima."
126
if isinstance(obj, types.FloatType):
128
elif isinstance(obj, types.IntType) or isinstance(obj, types.LongType):
130
elif isinstance(obj, MaskedArray) or isinstance(obj, ndarray):
132
if x in typecodes['Float']:
134
if x in typecodes['Integer']:
136
if x in typecodes['UnsignedInteger']:
139
raise TypeError, 'Unsuitable type for calculating maximum.'
141
def set_fill_value (a, fill_value):
142
"Set fill value of a if it is a masked array."
144
a.set_fill_value (fill_value)
147
"""Mask of values in a; could be nomask.
148
Returns nomask if a is not a masked array.
149
To get an array for sure use getmaskarray."""
150
if isinstance(a, MaskedArray):
155
def getmaskarray (a):
156
"""Mask of values in a; an array of zeros if mask is nomask
157
or not a masked array, and is a byte-sized integer.
158
Do not try to add up entries, for example.
162
return make_mask_none(shape(a))
167
"""Is m a legal mask? Does not check contents, only type.
170
return m.dtype.type is MaskType
171
except AttributeError:
174
def make_mask (m, copy=0, flag=0):
175
"""make_mask(m, copy=0, flag=0)
176
return m as a mask, creating a copy if necessary or requested.
177
Can accept any sequence of integers or nomask. Does not check
178
that contents must be 0s and 1s.
179
if flag, return nomask if m contains no true elements.
183
elif isinstance(m, ndarray):
184
if m.dtype.type is MaskType:
186
result = numeric.array(m, dtype=MaskType, copy=copy)
190
result = m.astype(MaskType)
192
result = filled(m, True).astype(MaskType)
194
if flag and not fromnumeric.sometrue(fromnumeric.ravel(result)):
199
def make_mask_none (s):
200
"Return a mask of all zeros of shape s."
201
result = numeric.zeros(s, dtype=MaskType)
205
def mask_or (m1, m2):
206
"""Logical or of the mask candidates m1 and m2, treating nomask as false.
207
Result may equal m1 or m2 if the other is nomask.
209
if m1 is nomask: return make_mask(m2)
210
if m2 is nomask: return make_mask(m1)
211
if m1 is m2 and is_mask(m1): return m1
212
return make_mask(umath.logical_or(m1, m2))
214
def filled (a, value = None):
215
"""a as a contiguous numeric array with any masked areas replaced by value
216
if value is None or the special element "masked", get_fill_value(a)
219
If a is already a contiguous numeric array, a itself is returned.
221
filled(a) can be used to be sure that the result is numeric when
222
passing an object a to other software ignorant of MA, in particular to
225
if isinstance(a, MaskedArray):
226
return a.filled(value)
227
elif isinstance(a, ndarray) and a.flags['CONTIGUOUS']:
229
elif isinstance(a, types.DictType):
230
return numeric.array(a, 'O')
232
return numeric.array(a)
234
def get_fill_value (a):
236
The fill value of a, if it has one; otherwise, the default fill value
240
result = a.fill_value()
242
result = default_fill_value(a)
245
def common_fill_value (a, b):
246
"The common fill_value of a and b, if there is one, or None"
247
t1 = get_fill_value(a)
248
t2 = get_fill_value(b)
249
if t1 == t2: return t1
252
# Domain functions return 1 where the argument(s) are not in the domain.
253
class domain_check_interval:
254
"domain_check_interval(a,b)(x) = true where x < a or y > b"
255
def __init__(self, y1, y2):
256
"domain_check_interval(a,b)(x) = true where x < a or y > b"
260
def __call__ (self, x):
261
"Execute the call behavior."
262
return umath.logical_or(umath.greater (x, self.y2),
263
umath.less(x, self.y1)
267
"domain_tan(eps) = true where abs(cos(x)) < eps)"
268
def __init__(self, eps):
269
"domain_tan(eps) = true where abs(cos(x)) < eps)"
272
def __call__ (self, x):
273
"Execute the call behavior."
274
return umath.less(umath.absolute(umath.cos(x)), self.eps)
276
class domain_greater:
277
"domain_greater(v)(x) = true where x <= v"
278
def __init__(self, critical_value):
279
"domain_greater(v)(x) = true where x <= v"
280
self.critical_value = critical_value
282
def __call__ (self, x):
283
"Execute the call behavior."
284
return umath.less_equal (x, self.critical_value)
286
class domain_greater_equal:
287
"domain_greater_equal(v)(x) = true where x < v"
288
def __init__(self, critical_value):
289
"domain_greater_equal(v)(x) = true where x < v"
290
self.critical_value = critical_value
292
def __call__ (self, x):
293
"Execute the call behavior."
294
return umath.less (x, self.critical_value)
296
class masked_unary_operation:
297
def __init__ (self, aufunc, fill=0, domain=None):
298
""" masked_unary_operation(aufunc, fill=0, domain=None)
299
aufunc(fill) must be defined
300
self(x) returns aufunc(x)
301
with masked values where domain(x) is true or getmask(x) is true.
306
self.__doc__ = getattr(aufunc, "__doc__", str(aufunc))
307
self.__name__ = getattr(aufunc, "__name__", str(aufunc))
308
ufunc_domain[aufunc] = domain
309
ufunc_fills[aufunc] = fill,
311
def __call__ (self, a, *args, **kwargs):
312
"Execute the call behavior."
313
# numeric tries to return scalars rather than arrays when given scalars.
315
d1 = filled(a, self.fill)
316
if self.domain is not None:
317
m = mask_or(m, self.domain(d1))
318
result = self.f(d1, *args, **kwargs)
319
return masked_array(result, m)
322
return "Masked version of " + str(self.f)
325
class domain_safe_divide:
326
def __init__ (self, tolerance=divide_tolerance):
327
self.tolerance = tolerance
328
def __call__ (self, a, b):
329
return umath.absolute(a) * self.tolerance >= umath.absolute(b)
331
class domained_binary_operation:
332
"""Binary operations that have a domain, like divide. These are complicated
333
so they are a separate class. They have no reduce, outer or accumulate.
335
def __init__ (self, abfunc, domain, fillx=0, filly=0):
336
"""abfunc(fillx, filly) must be defined.
337
abfunc(x, filly) = x for all x to enable reduce.
343
self.__doc__ = getattr(abfunc, "__doc__", str(abfunc))
344
self.__name__ = getattr(abfunc, "__name__", str(abfunc))
345
ufunc_domain[abfunc] = domain
346
ufunc_fills[abfunc] = fillx, filly
348
def __call__(self, a, b):
349
"Execute the call behavior."
352
d1 = filled(a, self.fillx)
353
d2 = filled(b, self.filly)
354
t = self.domain(d1, d2)
356
if fromnumeric.sometrue(t, None):
357
d2 = where(t, self.filly, d2)
360
result = self.f(d1, d2)
361
return masked_array(result, m)
364
return "Masked version of " + str(self.f)
366
class masked_binary_operation:
367
def __init__ (self, abfunc, fillx=0, filly=0):
368
"""abfunc(fillx, filly) must be defined.
369
abfunc(x, filly) = x for all x to enable reduce.
374
self.__doc__ = getattr(abfunc, "__doc__", str(abfunc))
375
ufunc_domain[abfunc] = None
376
ufunc_fills[abfunc] = fillx, filly
378
def __call__ (self, a, b, *args, **kwargs):
379
"Execute the call behavior."
380
m = mask_or(getmask(a), getmask(b))
381
d1 = filled(a, self.fillx)
382
d2 = filled(b, self.filly)
383
result = self.f(d1, d2, *args, **kwargs)
384
if isinstance(result, ndarray) \
386
and m.shape != result.shape:
387
m = mask_or(getmaskarray(a), getmaskarray(b))
388
return masked_array(result, m)
390
def reduce (self, target, axis=0, dtype=None):
391
"""Reduce target along the given axis with this function."""
393
t = filled(target, self.filly)
397
m = make_mask(m, copy=1)
400
t = self.f.reduce(t, axis)
402
t = masked_array (t, m)
403
# XXX: "or t.dtype" below is a workaround for what appears
404
# XXX: to be a bug in reduce.
405
t = self.f.reduce(filled(t, self.filly), axis,
406
dtype=dtype or t.dtype)
407
m = umath.logical_and.reduce(m, axis)
408
if isinstance(t, ndarray):
409
return masked_array(t, m, get_fill_value(target))
415
def outer (self, a, b):
416
"Return the function applied to the outer product of a and b."
419
if ma is nomask and mb is nomask:
424
m = logical_or.outer(ma, mb)
425
d = self.f.outer(filled(a, self.fillx), filled(b, self.filly))
426
return masked_array(d, m)
428
def accumulate (self, target, axis=0):
429
"""Accumulate target along axis after filling with y fill value."""
430
t = filled(target, self.filly)
431
return masked_array (self.f.accumulate (t, axis))
433
return "Masked version of " + str(self.f)
435
sqrt = masked_unary_operation(umath.sqrt, 0.0, domain_greater_equal(0.0))
436
log = masked_unary_operation(umath.log, 1.0, domain_greater(0.0))
437
log10 = masked_unary_operation(umath.log10, 1.0, domain_greater(0.0))
438
exp = masked_unary_operation(umath.exp)
439
conjugate = masked_unary_operation(umath.conjugate)
440
sin = masked_unary_operation(umath.sin)
441
cos = masked_unary_operation(umath.cos)
442
tan = masked_unary_operation(umath.tan, 0.0, domain_tan(1.e-35))
443
arcsin = masked_unary_operation(umath.arcsin, 0.0, domain_check_interval(-1.0, 1.0))
444
arccos = masked_unary_operation(umath.arccos, 0.0, domain_check_interval(-1.0, 1.0))
445
arctan = masked_unary_operation(umath.arctan)
446
# Missing from numeric
447
arcsinh = masked_unary_operation(umath.arcsinh)
448
arccosh = masked_unary_operation(umath.arccosh, 1.0, domain_greater_equal(1.0))
449
arctanh = masked_unary_operation(umath.arctanh, 0.0, domain_check_interval(-1.0+1e-15, 1.0-1e-15))
450
sinh = masked_unary_operation(umath.sinh)
451
cosh = masked_unary_operation(umath.cosh)
452
tanh = masked_unary_operation(umath.tanh)
453
absolute = masked_unary_operation(umath.absolute)
454
fabs = masked_unary_operation(umath.fabs)
455
negative = masked_unary_operation(umath.negative)
458
"""returns the indices of the elements of a which are not zero
461
return numeric.asarray(filled(a, 0).nonzero())
463
around = masked_unary_operation(fromnumeric.round_)
464
floor = masked_unary_operation(umath.floor)
465
ceil = masked_unary_operation(umath.ceil)
466
logical_not = masked_unary_operation(umath.logical_not)
468
add = masked_binary_operation(umath.add)
469
subtract = masked_binary_operation(umath.subtract)
470
subtract.reduce = None
471
multiply = masked_binary_operation(umath.multiply, 1, 1)
472
divide = domained_binary_operation(umath.divide, domain_safe_divide(), 0, 1)
473
true_divide = domained_binary_operation(umath.true_divide, domain_safe_divide(), 0, 1)
474
floor_divide = domained_binary_operation(umath.floor_divide, domain_safe_divide(), 0, 1)
475
remainder = domained_binary_operation(umath.remainder, domain_safe_divide(), 0, 1)
476
fmod = domained_binary_operation(umath.fmod, domain_safe_divide(), 0, 1)
477
hypot = masked_binary_operation(umath.hypot)
478
arctan2 = masked_binary_operation(umath.arctan2, 0.0, 1.0)
479
arctan2.reduce = None
480
equal = masked_binary_operation(umath.equal)
482
not_equal = masked_binary_operation(umath.not_equal)
483
not_equal.reduce = None
484
less_equal = masked_binary_operation(umath.less_equal)
485
less_equal.reduce = None
486
greater_equal = masked_binary_operation(umath.greater_equal)
487
greater_equal.reduce = None
488
less = masked_binary_operation(umath.less)
490
greater = masked_binary_operation(umath.greater)
491
greater.reduce = None
492
logical_and = masked_binary_operation(umath.logical_and)
493
alltrue = masked_binary_operation(umath.logical_and, 1, 1).reduce
494
logical_or = masked_binary_operation(umath.logical_or)
495
sometrue = logical_or.reduce
496
logical_xor = masked_binary_operation(umath.logical_xor)
497
bitwise_and = masked_binary_operation(umath.bitwise_and)
498
bitwise_or = masked_binary_operation(umath.bitwise_or)
499
bitwise_xor = masked_binary_operation(umath.bitwise_xor)
502
return fromnumeric.rank(filled(object))
505
return fromnumeric.shape(filled(object))
507
def size (object, axis=None):
508
return fromnumeric.size(filled(object), axis)
510
class MaskedArray (object):
511
"""Arrays with possibly masked values.
512
Masked values of 1 exclude the corresponding element from
516
x = array(data, dtype=None, copy=True, order=False,
517
mask = nomask, fill_value=None)
519
If copy=False, every effort is made not to copy the data:
520
If data is a MaskedArray, and argument mask=nomask,
521
then the candidate data is data.data and the
522
mask used is data.mask. If data is a numeric array,
523
it is used as the candidate raw data.
524
If dtype is not None and
525
is != data.dtype.char then a data copy is required.
526
Otherwise, the candidate is used.
528
If a data copy is required, raw data stored is the result of:
529
numeric.array(data, dtype=dtype.char, copy=copy)
531
If mask is nomask there are no masked values. Otherwise mask must
532
be convertible to an array of booleans with the same shape as x.
534
fill_value is used to fill in masked values when necessary,
535
such as when printing and in method/function filled().
536
The fill_value is not used for computation within this module.
538
__array_priority__ = 10.1
539
def __init__(self, data, dtype=None, copy=True, order=False,
540
mask=nomask, fill_value=None):
541
"""array(data, dtype=None, copy=True, order=False, mask=nomask, fill_value=None)
542
If data already a numeric array, its dtype becomes the default value of dtype.
547
tc = numeric.dtype(dtype)
548
need_data_copied = copy
549
if isinstance(data, MaskedArray):
554
need_data_copied = True
557
elif mask is not nomask: #attempting to change the mask
558
need_data_copied = True
560
elif isinstance(data, ndarray):
565
need_data_copied = True
567
need_data_copied = False #because I'll do it now
568
c = numeric.array(data, dtype=tc, copy=True, order=order)
573
self._data = numeric.array(c, dtype=tc, copy=True, order=order)
575
self._data = c.astype(tc)
581
self._shared_mask = 0
583
self._mask = make_mask (mask)
584
if self._mask is nomask:
585
self._shared_mask = 0
587
self._shared_mask = (self._mask is mask)
588
nm = size(self._mask)
589
nd = size(self._data)
592
self._mask = fromnumeric.resize(self._mask, self._data.shape)
593
self._shared_mask = 0
595
self._data = fromnumeric.resize(self._data, self._mask.shape)
596
self._data.shape = self._mask.shape
598
raise MAError, "Mask and data not compatible."
599
elif nm == 1 and shape(self._mask) != shape(self._data):
601
self._mask.shape = self._data.shape
603
self.set_fill_value(fill_value)
605
def __array__ (self, t=None, context=None):
606
"Special hook for numeric. Converts to numeric if possible."
607
if self._mask is not nomask:
608
if fromnumeric.ravel(self._mask).any():
610
warnings.warn("Cannot automatically convert masked array to "\
611
"numeric because data\n is masked in one or "\
615
# """Cannot automatically convert masked array to numeric because data
616
# is masked in one or more locations.
619
func, args, i = context
620
fills = ufunc_fills.get(func)
622
raise MAError, "%s not known to ma" % func
623
return self.filled(fills[i])
624
else: # Mask is all false
625
# Optimize to avoid future invocations of this section.
627
self._shared_mask = 0
629
return self._data.astype(t)
633
def __array_wrap__ (self, array, context=None):
634
"""Special hook for ufuncs.
636
Wraps the numpy array and sets the mask according to
640
return MaskedArray(array, copy=False, mask=nomask)
641
func, args = context[:2]
642
domain = ufunc_domain[func]
643
m = reduce(mask_or, [getmask(a) for a in args])
644
if domain is not None:
645
m = mask_or(m, domain(*[getattr(a, '_data', a)
650
except AttributeError:
654
m = reduce(mask_or, [getmaskarray(a) for a in args])
656
return MaskedArray(array, copy=False, mask=m)
658
def _get_shape(self):
659
"Return the current shape."
660
return self._data.shape
662
def _set_shape (self, newshape):
663
"Set the array's shape."
664
self._data.shape = newshape
665
if self._mask is not nomask:
666
self._mask = self._mask.copy()
667
self._mask.shape = newshape
670
"""Calculate the flat value.
672
if self._mask is nomask:
673
return masked_array(self._data.ravel(), mask=nomask,
674
fill_value = self.fill_value())
676
return masked_array(self._data.ravel(),
677
mask=self._mask.ravel(),
678
fill_value = self.fill_value())
680
def _set_flat (self, value):
686
"Get the real part of a complex array."
687
if self._mask is nomask:
688
return masked_array(self._data.real, mask=nomask,
689
fill_value = self.fill_value())
691
return masked_array(self._data.real, mask=self._mask,
692
fill_value = self.fill_value())
694
def _set_real (self, value):
699
def _get_imaginary(self):
700
"Get the imaginary part of a complex array."
701
if self._mask is nomask:
702
return masked_array(self._data.imag, mask=nomask,
703
fill_value = self.fill_value())
705
return masked_array(self._data.imag, mask=self._mask,
706
fill_value = self.fill_value())
708
def _set_imaginary (self, value):
709
"x.imaginary = value"
714
"""Calculate the str representation, using masked for fill if
715
it is enabled. Otherwise fill with fill value.
717
if masked_print_option.enabled():
718
f = masked_print_option
719
# XXX: Without the following special case masked
720
# XXX: would print as "[--]", not "--". Can we avoid
721
# XXX: checks for masked by choosing a different value
722
# XXX: for the masked singleton? 2005-01-05 -- sasha
726
if m is not nomask and m.shape == () and m:
728
# convert to object array to make filled work
729
self = self.astype(object)
731
f = self.fill_value()
736
"""Calculate the repr representation, using masked for fill if
737
it is enabled. Otherwise fill with fill value.
747
array(data = %(data)s,
751
without_mask = """array(
753
without_mask1 = """array(%(data)s)"""
756
if self._mask is nomask:
758
return without_mask1 % {'data':str(self.filled())}
759
return without_mask % {'data':str(self.filled())}
763
'data': str(self.filled()),
764
'mask': str(self._mask),
765
'fill': str(self.fill_value())
768
'data': str(self.filled()),
769
'mask': str(self._mask),
770
'fill': str(self.fill_value())
772
without_mask1 = """array(%(data)s)"""
773
if self._mask is nomask:
774
return without_mask % {'data':str(self.filled())}
777
'data': str(self.filled()),
778
'mask': str(self._mask),
779
'fill': str(self.fill_value())
783
"Convert self to float."
785
if self._mask is not nomask:
786
raise MAError, 'Cannot convert masked element to a Python float.'
787
return float(self.data.item())
790
"Convert self to int."
792
if self._mask is not nomask:
793
raise MAError, 'Cannot convert masked element to a Python int.'
794
return int(self.data.item())
796
def __getitem__(self, i):
797
"Get item described by i. Not a copy as in previous versions."
806
return masked_array(dout, fill_value=self._fill_value)
807
except AttributeError:
816
return masked_array(dout, mi, fill_value=self._fill_value)
819
# setitem and setslice notes
820
# note that if value is masked, it means to mask those locations.
821
# setting a value changes the mask to match the value in those locations.
823
def __setitem__(self, index, value):
824
"Set item described by index. If value is masked, mask those locations."
827
raise MAError, 'Cannot alter masked elements.'
829
if self._mask is nomask:
830
self._mask = make_mask_none(d.shape)
831
self._shared_mask = False
834
self._mask[index] = True
837
value = filled(value).astype(d.dtype)
840
if self._mask is not nomask:
842
self._mask[index] = False
844
if self._mask is nomask:
845
self._mask = make_mask_none(d.shape)
846
self._shared_mask = True
849
self._mask[index] = m
851
def __nonzero__(self):
852
"""returns true if any element is non-zero or masked
855
# XXX: This changes bool conversion logic from MA.
856
# XXX: In MA bool(a) == len(a) != 0, but in numpy
857
# XXX: scalars do not have len
860
return bool(m is not nomask and m.any()
861
or d is not nomask and d.any())
864
"""Return length of first dimension. This is weird but Python's
865
slicing behavior depends on it."""
866
return len(self._data)
868
def __and__(self, other):
870
return bitwise_and(self, other)
872
def __or__(self, other):
874
return bitwise_or(self, other)
876
def __xor__(self, other):
878
return bitwise_xor(self, other)
885
"Return absolute(self)"
886
return absolute(self)
889
"Return negative(self)"
890
return negative(self)
896
def __add__(self, other):
897
"Return add(self, other)"
898
return add(self, other)
902
def __mod__ (self, other):
903
"Return remainder(self, other)"
904
return remainder(self, other)
906
def __rmod__ (self, other):
907
"Return remainder(other, self)"
908
return remainder(other, self)
910
def __lshift__ (self, n):
911
return left_shift(self, n)
913
def __rshift__ (self, n):
914
return right_shift(self, n)
916
def __sub__(self, other):
917
"Return subtract(self, other)"
918
return subtract(self, other)
920
def __rsub__(self, other):
921
"Return subtract(other, self)"
922
return subtract(other, self)
924
def __mul__(self, other):
925
"Return multiply(self, other)"
926
return multiply(self, other)
930
def __div__(self, other):
931
"Return divide(self, other)"
932
return divide(self, other)
934
def __rdiv__(self, other):
935
"Return divide(other, self)"
936
return divide(other, self)
938
def __truediv__(self, other):
939
"Return divide(self, other)"
940
return true_divide(self, other)
942
def __rtruediv__(self, other):
943
"Return divide(other, self)"
944
return true_divide(other, self)
946
def __floordiv__(self, other):
947
"Return divide(self, other)"
948
return floor_divide(self, other)
950
def __rfloordiv__(self, other):
951
"Return divide(other, self)"
952
return floor_divide(other, self)
954
def __pow__(self, other, third=None):
955
"Return power(self, other, third)"
956
return power(self, other, third)
962
def __iadd__(self, other):
963
"Add other to self in place."
964
t = self._data.dtype.char
969
elif t in typecodes['Integer']:
970
if t1 in typecodes['Integer']:
973
raise TypeError, 'Incorrect type for in-place operation.'
974
elif t in typecodes['Float']:
975
if t1 in typecodes['Integer']:
977
elif t1 in typecodes['Float']:
980
raise TypeError, 'Incorrect type for in-place operation.'
981
elif t in typecodes['Complex']:
982
if t1 in typecodes['Integer']:
984
elif t1 in typecodes['Float']:
986
elif t1 in typecodes['Complex']:
989
raise TypeError, 'Incorrect type for in-place operation.'
991
raise TypeError, 'Incorrect type for in-place operation.'
993
if self._mask is nomask:
997
self._shared_mask = m is not nomask
999
result = add(self, masked_array(f, mask=getmask(other)))
1000
self._data = result.data
1001
self._mask = result.mask
1002
self._shared_mask = 1
1005
def __imul__(self, other):
1006
"Add other to self in place."
1007
t = self._data.dtype.char
1008
f = filled(other, 0)
1012
elif t in typecodes['Integer']:
1013
if t1 in typecodes['Integer']:
1016
raise TypeError, 'Incorrect type for in-place operation.'
1017
elif t in typecodes['Float']:
1018
if t1 in typecodes['Integer']:
1020
elif t1 in typecodes['Float']:
1023
raise TypeError, 'Incorrect type for in-place operation.'
1024
elif t in typecodes['Complex']:
1025
if t1 in typecodes['Integer']:
1027
elif t1 in typecodes['Float']:
1029
elif t1 in typecodes['Complex']:
1032
raise TypeError, 'Incorrect type for in-place operation.'
1034
raise TypeError, 'Incorrect type for in-place operation.'
1036
if self._mask is nomask:
1040
self._shared_mask = m is not nomask
1042
result = multiply(self, masked_array(f, mask=getmask(other)))
1043
self._data = result.data
1044
self._mask = result.mask
1045
self._shared_mask = 1
1048
def __isub__(self, other):
1049
"Subtract other from self in place."
1050
t = self._data.dtype.char
1051
f = filled(other, 0)
1055
elif t in typecodes['Integer']:
1056
if t1 in typecodes['Integer']:
1059
raise TypeError, 'Incorrect type for in-place operation.'
1060
elif t in typecodes['Float']:
1061
if t1 in typecodes['Integer']:
1063
elif t1 in typecodes['Float']:
1066
raise TypeError, 'Incorrect type for in-place operation.'
1067
elif t in typecodes['Complex']:
1068
if t1 in typecodes['Integer']:
1070
elif t1 in typecodes['Float']:
1072
elif t1 in typecodes['Complex']:
1075
raise TypeError, 'Incorrect type for in-place operation.'
1077
raise TypeError, 'Incorrect type for in-place operation.'
1079
if self._mask is nomask:
1083
self._shared_mask = m is not nomask
1085
result = subtract(self, masked_array(f, mask=getmask(other)))
1086
self._data = result.data
1087
self._mask = result.mask
1088
self._shared_mask = 1
1093
def __idiv__(self, other):
1094
"Divide self by other in place."
1095
t = self._data.dtype.char
1096
f = filled(other, 0)
1100
elif t in typecodes['Integer']:
1101
if t1 in typecodes['Integer']:
1104
raise TypeError, 'Incorrect type for in-place operation.'
1105
elif t in typecodes['Float']:
1106
if t1 in typecodes['Integer']:
1108
elif t1 in typecodes['Float']:
1111
raise TypeError, 'Incorrect type for in-place operation.'
1112
elif t in typecodes['Complex']:
1113
if t1 in typecodes['Integer']:
1115
elif t1 in typecodes['Float']:
1117
elif t1 in typecodes['Complex']:
1120
raise TypeError, 'Incorrect type for in-place operation.'
1122
raise TypeError, 'Incorrect type for in-place operation.'
1124
result = divide(self, masked_array(f, mask=mo))
1125
self._data = result.data
1126
dm = result.raw_mask()
1127
if dm is not self._mask:
1129
self._shared_mask = 1
1132
def __eq__(self, other):
1133
return equal(self,other)
1135
def __ne__(self, other):
1136
return not_equal(self,other)
1138
def __lt__(self, other):
1139
return less(self,other)
1141
def __le__(self, other):
1142
return less_equal(self,other)
1144
def __gt__(self, other):
1145
return greater(self,other)
1147
def __ge__(self, other):
1148
return greater_equal(self,other)
1150
def astype (self, tc):
1151
"return self as array of given type."
1152
d = self._data.astype(tc)
1153
return array(d, mask=self._mask)
1155
def byte_swapped(self):
1156
"""Returns the raw data field, byte_swapped. Included for consistency
1157
with numeric but doesn't make sense in this context.
1159
return self._data.byte_swapped()
1161
def compressed (self):
1162
"A 1-D array of all the non-masked data."
1163
d = fromnumeric.ravel(self._data)
1164
if self._mask is nomask:
1167
m = 1 - fromnumeric.ravel(self._mask)
1168
c = fromnumeric.compress(m, d)
1169
return array(c, copy=0)
1171
def count (self, axis = None):
1172
"Count of the non-masked elements in a, or along a certain axis."
1174
s = self._data.shape
1182
return reduce(lambda x, y:x*y, s)
1189
w = fromnumeric.ravel(m).astype(int)
1194
n2 = umath.add.reduce(w)
1198
n2 = sum(m.astype(int), axis)
1201
def dot (self, other):
1202
"s.dot(other) = innerproduct(s, other)"
1203
return innerproduct(self, other)
1205
def fill_value(self):
1206
"Get the current fill value."
1207
return self._fill_value
1209
def filled (self, fill_value=None):
1210
"""A numeric array with masked values filled. If fill_value is None,
1211
use self.fill_value().
1213
If mask is nomask, copy data only if not contiguous.
1214
Result is always a contiguous, numeric array.
1215
# Is contiguous really necessary now?
1220
if d.flags['CONTIGUOUS']:
1225
if fill_value is None:
1226
value = self._fill_value
1231
result = numeric.array(value)
1234
result = numeric.array(d, dtype=d.dtype, copy=1)
1236
except (TypeError, AttributeError):
1237
#ok, can't put that value in here
1238
value = numeric.array(value, dtype=object)
1239
d = d.astype(object)
1240
result = fromnumeric.choose(m, (d, value))
1244
"""Return the ids of the data and mask areas"""
1245
return (id(self._data), id(self._mask))
1247
def iscontiguous (self):
1248
"Is the data contiguous?"
1249
return self._data.flags['CONTIGUOUS']
1252
"Item size of each data item."
1253
return self._data.itemsize
1256
def outer(self, other):
1257
"s.outer(other) = outerproduct(s, other)"
1258
return outerproduct(self, other)
1260
def put (self, values):
1261
"""Set the non-masked entries of self to filled(values).
1264
iota = numeric.arange(self.size)
1266
if self._mask is nomask:
1269
ind = fromnumeric.compress(1 - self._mask, iota)
1270
d[ind] = filled(values).astype(d.dtype)
1272
def putmask (self, values):
1273
"""Set the masked entries of self to filled(values).
1274
Mask changed to nomask.
1277
if self._mask is not nomask:
1278
d[self._mask] = filled(values).astype(d.dtype)
1279
self._shared_mask = 0
1283
"""Return a 1-D view of self."""
1284
if self._mask is nomask:
1285
return masked_array(self._data.ravel())
1287
return masked_array(self._data.ravel(), self._mask.ravel())
1289
def raw_data (self):
1290
""" Obsolete; use data property instead.
1291
The raw data; portions may be meaningless.
1292
May be noncontiguous. Expert use only."""
1294
data = property(fget=raw_data,
1295
doc="The data, but values at masked locations are meaningless.")
1297
def raw_mask (self):
1298
""" Obsolete; use mask property instead.
1299
May be noncontiguous. Expert use only.
1302
mask = property(fget=raw_mask,
1303
doc="The mask, may be nomask. Values where mask true are meaningless.")
1305
def reshape (self, *s):
1306
"""This array reshaped to shape s"""
1307
d = self._data.reshape(*s)
1308
if self._mask is nomask:
1309
return masked_array(d)
1311
m = self._mask.reshape(*s)
1312
return masked_array(d, m)
1314
def set_fill_value (self, v=None):
1315
"Set the fill value to v. Omit v to restore default."
1317
v = default_fill_value (self.raw_data())
1318
self._fill_value = v
1320
def _get_ndim(self):
1321
return self._data.ndim
1322
ndim = property(_get_ndim, doc=numeric.ndarray.ndim.__doc__)
1324
def _get_size (self):
1325
return self._data.size
1326
size = property(fget=_get_size, doc="Number of elements in the array.")
1327
## CHECK THIS: signature of numeric.array.size?
1329
def _get_dtype(self):
1330
return self._data.dtype
1331
dtype = property(fget=_get_dtype, doc="type of the array elements.")
1333
def item(self, *args):
1334
"Return Python scalar if possible"
1335
if self._mask is not nomask:
1336
m = self._mask.item(*args)
1342
return self._data.item(*args)
1344
def itemset(self, *args):
1345
"Set Python scalar into array"
1350
def tolist(self, fill_value=None):
1352
return self.filled(fill_value).tolist()
1354
def tostring(self, fill_value=None):
1356
return self.filled(fill_value).tostring()
1359
"Replace the mask by nomask if possible."
1360
if self._mask is nomask: return
1361
m = make_mask(self._mask, flag=1)
1364
self._shared_mask = 0
1366
def unshare_mask (self):
1367
"If currently sharing mask, make a copy."
1368
if self._shared_mask:
1369
self._mask = make_mask (self._mask, copy=1, flag=0)
1370
self._shared_mask = 0
1372
def _get_ctypes(self):
1373
return self._data.ctypes
1378
return self.transpose()
1380
shape = property(_get_shape, _set_shape,
1381
doc = 'tuple giving the shape of the array')
1383
flat = property(_get_flat, _set_flat,
1384
doc = 'Access array in flat form.')
1386
real = property(_get_real, _set_real,
1387
doc = 'Access the real part of the array')
1389
imaginary = property(_get_imaginary, _set_imaginary,
1390
doc = 'Access the imaginary part of the array')
1394
ctypes = property(_get_ctypes, None, doc="ctypes")
1396
T = property(_get_T, None, doc="get transpose")
1398
#end class MaskedArray
1402
def isMaskedArray (x):
1403
"Is x a masked array, that is, an instance of MaskedArray?"
1404
return isinstance(x, MaskedArray)
1406
isarray = isMaskedArray
1407
isMA = isMaskedArray #backward compatibility
1409
def allclose (a, b, fill_value=1, rtol=1.e-5, atol=1.e-8):
1410
""" Returns true if all components of a and b are equal
1411
subject to given tolerances.
1412
If fill_value is 1, masked values considered equal.
1413
If fill_value is 0, masked values considered unequal.
1414
The relative error rtol should be positive and << 1.0
1415
The absolute error atol comes into play for those elements
1416
of b that are very small or zero; it says how small a must be also.
1418
m = mask_or(getmask(a), getmask(b))
1421
x = filled(array(d1, copy=0, mask=m), fill_value).astype(float)
1422
y = filled(array(d2, copy=0, mask=m), 1).astype(float)
1423
d = umath.less_equal(umath.absolute(x-y), atol + rtol * umath.absolute(y))
1424
return fromnumeric.alltrue(fromnumeric.ravel(d))
1426
def allequal (a, b, fill_value=1):
1428
True if all entries of a and b are equal, using
1429
fill_value as a truth value where either or both are masked.
1431
m = mask_or(getmask(a), getmask(b))
1435
d = umath.equal(x, y)
1436
return fromnumeric.alltrue(fromnumeric.ravel(d))
1440
d = umath.equal(x, y)
1441
dm = array(d, mask=m, copy=0)
1442
return fromnumeric.alltrue(fromnumeric.ravel(filled(dm, 1)))
1446
def masked_values (data, value, rtol=1.e-5, atol=1.e-8, copy=1):
1448
masked_values(data, value, rtol=1.e-5, atol=1.e-8)
1449
Create a masked array; mask is nomask if possible.
1450
If copy==0, and otherwise possible, result
1451
may share data values with original array.
1452
Let d = filled(data, value). Returns d
1453
masked where abs(data-value)<= atol + rtol * abs(value)
1454
if d is of a floating point type. Otherwise returns
1455
masked_object(d, value, copy)
1457
abs = umath.absolute
1458
d = filled(data, value)
1459
if issubclass(d.dtype.type, numeric.floating):
1460
m = umath.less_equal(abs(d-value), atol+rtol*abs(value))
1461
m = make_mask(m, flag=1)
1462
return array(d, mask = m, copy=copy,
1465
return masked_object(d, value, copy=copy)
1467
def masked_object (data, value, copy=1):
1468
"Create array masked where exactly data equal to value"
1469
d = filled(data, value)
1470
dm = make_mask(umath.equal(d, value), flag=1)
1471
return array(d, mask=dm, copy=copy, fill_value=value)
1473
def arange(start, stop=None, step=1, dtype=None):
1474
"""Just like range() except it returns a array whose type can be specified
1475
by the keyword argument dtype.
1477
return array(numeric.arange(start, stop, step, dtype))
1481
def fromstring (s, t):
1482
"Construct a masked array from a string. Result will have no mask."
1483
return masked_array(numeric.fromstring(s, t))
1485
def left_shift (a, n):
1489
d = umath.left_shift(filled(a), n)
1490
return masked_array(d)
1492
d = umath.left_shift(filled(a, 0), n)
1493
return masked_array(d, m)
1495
def right_shift (a, n):
1496
"Right shift n bits"
1499
d = umath.right_shift(filled(a), n)
1500
return masked_array(d)
1502
d = umath.right_shift(filled(a, 0), n)
1503
return masked_array(d, m)
1505
def resize (a, new_shape):
1506
"""resize(a, new_shape) returns a new array with the specified shape.
1507
The original array's total size can be any size."""
1510
m = fromnumeric.resize(m, new_shape)
1511
result = array(fromnumeric.resize(filled(a), new_shape), mask=m)
1512
result.set_fill_value(get_fill_value(a))
1515
def new_repeat(a, repeats, axis=None):
1516
"""repeat elements of a repeats times along axis
1517
repeats is a sequence of length a.shape[axis]
1518
telling how many times to repeat each element.
1521
if isinstance(repeats, types.IntType):
1525
num = af.shape[axis]
1526
repeats = tuple([repeats]*num)
1530
m = fromnumeric.repeat(m, repeats, axis)
1531
d = fromnumeric.repeat(af, repeats, axis)
1532
result = masked_array(d, m)
1533
result.set_fill_value(get_fill_value(a))
1539
"""identity(n) returns the identity matrix of shape n x n.
1541
return array(numeric.identity(n))
1543
def indices (dimensions, dtype=None):
1544
"""indices(dimensions,dtype=None) returns an array representing a grid
1545
of indices with row-only, and column-only variation.
1547
return array(numeric.indices(dimensions, dtype))
1549
def zeros (shape, dtype=float):
1550
"""zeros(n, dtype=float) =
1551
an array of all zeros of the given length or shape."""
1552
return array(numeric.zeros(shape, dtype))
1554
def ones (shape, dtype=float):
1555
"""ones(n, dtype=float) =
1556
an array of all ones of the given length or shape."""
1557
return array(numeric.ones(shape, dtype))
1559
def count (a, axis = None):
1560
"Count of the non-masked elements in a, or along a certain axis."
1562
return a.count(axis)
1564
def power (a, b, third=None):
1566
if third is not None:
1567
raise MAError, "3-argument power not supported."
1573
if fb.dtype.char in typecodes["Integer"]:
1574
return masked_array(umath.power(fa, fb), m)
1575
md = make_mask(umath.less(fa, 0), flag=1)
1578
return masked_array(umath.power(fa, fb))
1580
fa = numeric.where(m, 1, fa)
1581
return masked_array(umath.power(fa, fb), m)
1583
def masked_array (a, mask=nomask, fill_value=None):
1584
"""masked_array(a, mask=nomask) =
1585
array(a, mask=mask, copy=0, fill_value=fill_value)
1587
return array(a, mask=mask, copy=0, fill_value=fill_value)
1589
def sum (target, axis=None, dtype=None):
1591
target = ravel(target)
1593
return add.reduce(target, axis, dtype)
1595
def product (target, axis=None, dtype=None):
1597
target = ravel(target)
1599
return multiply.reduce(target, axis, dtype)
1601
def new_average (a, axis=None, weights=None, returned = 0):
1602
"""average(a, axis=None, weights=None)
1603
Computes average along indicated axis.
1604
If axis is None, average over the entire array
1605
Inputs can be integer or floating types; result is of type float.
1607
If weights are given, result is sum(a*weights,axis=0)/(sum(weights,axis=0)*1.0)
1608
weights must have a's shape or be the 1-d with length the size
1609
of a in the given axis.
1611
If returned, return a tuple: the result and the sum of the weights
1612
or count of values. Results will have the same shape.
1614
masked values in the weights will be set to 0.0
1624
n = add.reduce(a.raw_data().ravel())
1625
d = reduce(lambda x, y: x * y, ash, 1.0)
1627
w = filled(weights, 0.0).ravel()
1628
n = umath.add.reduce(a.raw_data().ravel() * w)
1629
d = umath.add.reduce(w)
1633
n = add.reduce(a.ravel())
1634
w = fromnumeric.choose(mask, (1.0, 0.0)).ravel()
1635
d = umath.add.reduce(w)
1638
w = array(filled(weights, 0.0), float, mask=mask).ravel()
1639
n = add.reduce(a.ravel() * w)
1646
n = umath.add.reduce(a.raw_data(), axis)
1648
w = filled(weights, 0.0)
1653
w = numeric.array(w, float, copy=0)
1654
n = add.reduce(a*w, axis)
1655
d = add.reduce(w, axis)
1657
elif wsh == (ash[axis],):
1658
r = [newaxis]*len(ash)
1659
r[axis] = slice(None, None, 1)
1660
w = eval ("w["+ repr(tuple(r)) + "] * ones(ash, float)")
1661
n = add.reduce(a*w, axis)
1662
d = add.reduce(w, axis)
1665
raise ValueError, 'average: weights wrong shape.'
1668
n = add.reduce(a, axis)
1669
w = numeric.choose(mask, (1.0, 0.0))
1670
d = umath.add.reduce(w, axis)
1673
w = filled(weights, 0.0)
1678
w = array(w, float, mask=mask, copy=0)
1679
n = add.reduce(a*w, axis)
1680
d = add.reduce(w, axis)
1681
elif wsh == (ash[axis],):
1682
r = [newaxis]*len(ash)
1683
r[axis] = slice(None, None, 1)
1684
w = eval ("w["+ repr(tuple(r)) + "] * masked_array(ones(ash, float), mask)")
1685
n = add.reduce(a*w, axis)
1686
d = add.reduce(w, axis)
1688
raise ValueError, 'average: weights wrong shape.'
1690
#print n, d, repr(mask), repr(weights)
1691
if n is masked or d is masked: return masked
1692
result = divide (n, d)
1695
if isinstance(result, MaskedArray):
1698
if not isinstance(d, MaskedArray):
1700
if not d.shape == result.shape:
1701
d = ones(result.shape, float) * d
1708
def where (condition, x, y):
1709
"""where(condition, x, y) is x where condition is nonzero, y otherwise.
1710
condition must be convertible to an integer array.
1711
Answer is always the shape of condition.
1712
The type depends on x and y. It is integer if both x and y are
1715
fc = filled(not_equal(condition, 0), 0)
1720
d = numeric.choose(fc, (yv, xv))
1721
md = numeric.choose(fc, (ym, xm))
1722
m = getmask(condition)
1723
m = make_mask(mask_or(m, md), copy=0, flag=1)
1724
return masked_array(d, m)
1726
def choose (indices, t, out=None, mode='raise'):
1727
"Returns array shaped like indices with elements chosen from t"
1729
if x is masked: return 1
1732
if x is masked: return 1
1734
if m is nomask: return 0
1736
c = filled(indices, 0)
1737
masks = [nmask(x) for x in t]
1738
a = [fmask(x) for x in t]
1739
d = numeric.choose(c, a)
1740
m = numeric.choose(c, masks)
1741
m = make_mask(mask_or(m, getmask(indices)), copy=0, flag=1)
1742
return masked_array(d, m)
1744
def masked_where(condition, x, copy=1):
1745
"""Return x as an array masked where condition is true.
1746
Also masked where x or condition masked.
1748
cm = filled(condition,1)
1749
m = mask_or(getmask(x), cm)
1750
return array(filled(x), copy=copy, mask=m)
1752
def masked_greater(x, value, copy=1):
1753
"masked_greater(x, value) = x masked where x > value"
1754
return masked_where(greater(x, value), x, copy)
1756
def masked_greater_equal(x, value, copy=1):
1757
"masked_greater_equal(x, value) = x masked where x >= value"
1758
return masked_where(greater_equal(x, value), x, copy)
1760
def masked_less(x, value, copy=1):
1761
"masked_less(x, value) = x masked where x < value"
1762
return masked_where(less(x, value), x, copy)
1764
def masked_less_equal(x, value, copy=1):
1765
"masked_less_equal(x, value) = x masked where x <= value"
1766
return masked_where(less_equal(x, value), x, copy)
1768
def masked_not_equal(x, value, copy=1):
1769
"masked_not_equal(x, value) = x masked where x != value"
1771
c = umath.not_equal(d, value)
1772
m = mask_or(c, getmask(x))
1773
return array(d, mask=m, copy=copy)
1775
def masked_equal(x, value, copy=1):
1776
"""masked_equal(x, value) = x masked where x == value
1777
For floating point consider masked_values(x, value) instead.
1780
c = umath.equal(d, value)
1781
m = mask_or(c, getmask(x))
1782
return array(d, mask=m, copy=copy)
1784
def masked_inside(x, v1, v2, copy=1):
1785
"""x with mask of all values of x that are inside [v1,v2]
1786
v1 and v2 can be given in either order.
1793
c = umath.logical_and(umath.less_equal(d, v2), umath.greater_equal(d, v1))
1794
m = mask_or(c, getmask(x))
1795
return array(d, mask = m, copy=copy)
1797
def masked_outside(x, v1, v2, copy=1):
1798
"""x with mask of all values of x that are outside [v1,v2]
1799
v1 and v2 can be given in either order.
1806
c = umath.logical_or(umath.less(d, v1), umath.greater(d, v2))
1807
m = mask_or(c, getmask(x))
1808
return array(d, mask = m, copy=copy)
1810
def reshape (a, *newshape):
1811
"Copy of a with a new shape."
1813
d = filled(a).reshape(*newshape)
1815
return masked_array(d)
1817
return masked_array(d, mask=numeric.reshape(m, *newshape))
1820
"a as one-dimensional, may share data and mask"
1822
d = fromnumeric.ravel(filled(a))
1824
return masked_array(d)
1826
return masked_array(d, mask=numeric.ravel(m))
1828
def concatenate (arrays, axis=0):
1829
"Concatenate the arrays along the given axis"
1833
d = numeric.concatenate(d, axis)
1835
if getmask(x) is not nomask: break
1837
return masked_array(d)
1840
dm.append(getmaskarray(x))
1841
dm = numeric.concatenate(dm, axis)
1842
return masked_array(d, mask=dm)
1844
def swapaxes (a, axis1, axis2):
1846
d = masked_array(a).data
1848
return masked_array(data=numeric.swapaxes(d, axis1, axis2))
1850
return masked_array(data=numeric.swapaxes(d, axis1, axis2),
1851
mask=numeric.swapaxes(m, axis1, axis2),)
1854
def new_take (a, indices, axis=None, out=None, mode='raise'):
1855
"returns selection of items from a."
1857
# d = masked_array(a).raw_data()
1858
d = masked_array(a).data
1860
return masked_array(numeric.take(d, indices, axis))
1862
return masked_array(numeric.take(d, indices, axis),
1863
mask = numeric.take(m, indices, axis))
1865
def transpose(a, axes=None):
1866
"reorder dimensions per tuple axes"
1870
return masked_array(numeric.transpose(d, axes))
1872
return masked_array(numeric.transpose(d, axes),
1873
mask = numeric.transpose(m, axes))
1876
def put(a, indices, values, mode='raise'):
1877
"""sets storage-indexed locations to corresponding values.
1879
Values and indices are filled if necessary.
1883
ind = filled(indices)
1885
numeric.put (d, ind, v)
1889
numeric.put(a.raw_mask(), ind, 0)
1891
def putmask(a, mask, values):
1892
"putmask(a, mask, values) sets a where mask is true."
1895
numeric.putmask(a.raw_data(), mask, values)
1897
if m is nomask: return
1899
numeric.putmask(a.raw_mask(), mask, 0)
1902
"""inner(a,b) returns the dot product of two arrays, which has
1903
shape a.shape[:-1] + b.shape[:-1] with elements computed by summing the
1904
product of the elements from the last dimensions of a and b.
1905
Masked elements are replace by zeros.
1909
if len(fa.shape) == 0: fa.shape = (1,)
1910
if len(fb.shape) == 0: fb.shape = (1,)
1911
return masked_array(numeric.inner(fa, fb))
1913
innerproduct = inner
1916
"""outer(a,b) = {a[i]*b[j]}, has shape (len(a),len(b))"""
1917
fa = filled(a, 0).ravel()
1918
fb = filled(b, 0).ravel()
1919
d = numeric.outer(fa, fb)
1922
if ma is nomask and mb is nomask:
1923
return masked_array(d)
1924
ma = getmaskarray(a)
1925
mb = getmaskarray(b)
1926
m = make_mask(1-numeric.outer(1-ma, 1-mb), copy=0)
1927
return masked_array(d, m)
1929
outerproduct = outer
1932
"""dot(a,b) returns matrix-multiplication between a and b. The product-sum
1933
is over the last dimension of a and the second-to-last dimension of b.
1934
Masked values are replaced by zeros. See also innerproduct.
1936
return innerproduct(filled(a, 0), numeric.swapaxes(filled(b, 0), -1, -2))
1938
def compress(condition, x, dimension=-1, out=None):
1939
"""Select those parts of x for which condition is true.
1940
Masked values in condition are considered false.
1942
c = filled(condition, 0)
1945
m = numeric.compress(c, m, dimension)
1946
d = numeric.compress(c, filled(x), dimension)
1947
return masked_array(d, m)
1949
class _minimum_operation:
1950
"Object to calculate minima"
1951
def __init__ (self):
1952
"""minimum(a, b) or minimum(a)
1953
In one argument case returns the scalar minimum.
1957
def __call__ (self, a, b=None):
1958
"Execute the call behavior."
1962
d = amin(filled(a).ravel())
1968
return amin(ac.raw_data())
1970
return where(less(a, b), a, b)
1972
def reduce (self, target, axis=0):
1973
"""Reduce target along the given axis."""
1977
return masked_array (umath.minimum.reduce (t, axis))
1979
t = umath.minimum.reduce(filled(target, minimum_fill_value(target)), axis)
1980
m = umath.logical_and.reduce(m, axis)
1981
return masked_array(t, m, get_fill_value(target))
1983
def outer (self, a, b):
1984
"Return the function applied to the outer product of a and b."
1987
if ma is nomask and mb is nomask:
1990
ma = getmaskarray(a)
1991
mb = getmaskarray(b)
1992
m = logical_or.outer(ma, mb)
1993
d = umath.minimum.outer(filled(a), filled(b))
1994
return masked_array(d, m)
1996
minimum = _minimum_operation ()
1998
class _maximum_operation:
1999
"Object to calculate maxima"
2000
def __init__ (self):
2001
"""maximum(a, b) or maximum(a)
2002
In one argument case returns the scalar maximum.
2006
def __call__ (self, a, b=None):
2007
"Execute the call behavior."
2011
d = amax(filled(a).ravel())
2017
return amax(ac.raw_data())
2019
return where(greater(a, b), a, b)
2021
def reduce (self, target, axis=0):
2022
"""Reduce target along the given axis."""
2026
return masked_array (umath.maximum.reduce (t, axis))
2028
t = umath.maximum.reduce(filled(target, maximum_fill_value(target)), axis)
2029
m = umath.logical_and.reduce(m, axis)
2030
return masked_array(t, m, get_fill_value(target))
2032
def outer (self, a, b):
2033
"Return the function applied to the outer product of a and b."
2036
if ma is nomask and mb is nomask:
2039
ma = getmaskarray(a)
2040
mb = getmaskarray(b)
2041
m = logical_or.outer(ma, mb)
2042
d = umath.maximum.outer(filled(a), filled(b))
2043
return masked_array(d, m)
2045
maximum = _maximum_operation ()
2047
def sort (x, axis = -1, fill_value=None):
2048
"""If x does not have a mask, return a masked array formed from the
2049
result of numeric.sort(x, axis).
2050
Otherwise, fill x with fill_value. Sort it.
2051
Set a mask where the result is equal to fill_value.
2052
Note that this may have unintended consequences if the data contains the
2053
fill value at a non-masked site.
2055
If fill_value is not given the default fill value for x's type will be
2058
if fill_value is None:
2059
fill_value = default_fill_value (x)
2060
d = filled(x, fill_value)
2061
s = fromnumeric.sort(d, axis)
2062
if getmask(x) is nomask:
2063
return masked_array(s)
2064
return masked_values(s, fill_value, copy=0)
2066
def diagonal(a, k = 0, axis1=0, axis2=1):
2067
"""diagonal(a,k=0,axis1=0, axis2=1) = the k'th diagonal of a"""
2068
d = fromnumeric.diagonal(filled(a), k, axis1, axis2)
2071
return masked_array(d, m)
2073
return masked_array(d, fromnumeric.diagonal(m, k, axis1, axis2))
2075
def trace (a, offset=0, axis1=0, axis2=1, dtype=None, out=None):
2076
"""trace(a,offset=0, axis1=0, axis2=1) returns the sum along diagonals
2077
(defined by the last two dimenions) of the array.
2079
return diagonal(a, offset, axis1, axis2).sum(dtype=dtype)
2081
def argsort (x, axis = -1, out=None, fill_value=None):
2082
"""Treating masked values as if they have the value fill_value,
2083
return sort indices for sorting along given axis.
2084
if fill_value is None, use get_fill_value(x)
2085
Returns a numpy array.
2087
d = filled(x, fill_value)
2088
return fromnumeric.argsort(d, axis)
2090
def argmin (x, axis = -1, out=None, fill_value=None):
2091
"""Treating masked values as if they have the value fill_value,
2092
return indices for minimum values along given axis.
2093
if fill_value is None, use get_fill_value(x).
2094
Returns a numpy array if x has more than one dimension.
2095
Otherwise, returns a scalar index.
2097
d = filled(x, fill_value)
2098
return fromnumeric.argmin(d, axis)
2100
def argmax (x, axis = -1, out=None, fill_value=None):
2101
"""Treating masked values as if they have the value fill_value,
2102
return sort indices for maximum along given axis.
2103
if fill_value is None, use -get_fill_value(x) if it exists.
2104
Returns a numpy array if x has more than one dimension.
2105
Otherwise, returns a scalar index.
2107
if fill_value is None:
2108
fill_value = default_fill_value (x)
2110
fill_value = - fill_value
2113
d = filled(x, fill_value)
2114
return fromnumeric.argmax(d, axis)
2116
def fromfunction (f, s):
2117
"""apply f to s to create array as in umath."""
2118
return masked_array(numeric.fromfunction(f, s))
2120
def asarray(data, dtype=None):
2121
"""asarray(data, dtype) = array(data, dtype, copy=0)
2123
if isinstance(data, MaskedArray) and \
2124
(dtype is None or dtype == data.dtype):
2126
return array(data, dtype=dtype, copy=0)
2128
# Add methods to support ndarray interface
2129
# XXX: I is better to to change the masked_*_operation adaptors
2130
# XXX: to wrap ndarray methods directly to create ma.array methods.
2131
from types import MethodType
2133
return MethodType(f, None, array)
2134
def not_implemented(*args, **kwds):
2135
raise NotImplementedError, "not yet implemented for numpy.ma arrays"
2136
array.all = _m(alltrue)
2137
array.any = _m(sometrue)
2138
array.argmax = _m(argmax)
2139
array.argmin = _m(argmin)
2140
array.argsort = _m(argsort)
2141
array.base = property(_m(not_implemented))
2142
array.byteswap = _m(not_implemented)
2144
def _choose(self, *args, **kwds):
2145
return choose(self, args)
2146
array.choose = _m(_choose)
2149
def _clip(self,a_min,a_max,out=None):
2150
return MaskedArray(data = self.data.clip(asarray(a_min).data,
2151
asarray(a_max).data),
2152
mask = mask_or(self.mask,
2153
mask_or(getmask(a_min),getmask(a_max))))
2154
array.clip = _m(_clip)
2156
def _compress(self, cond, axis=None, out=None):
2157
return compress(cond, self, axis)
2158
array.compress = _m(_compress)
2161
array.conj = array.conjugate = _m(conjugate)
2162
array.copy = _m(not_implemented)
2164
def _cumprod(self, axis=None, dtype=None, out=None):
2167
m = umath.logical_or.accumulate(self.mask, axis)
2168
return MaskedArray(data = self.filled(1).cumprod(axis, dtype), mask=m)
2169
array.cumprod = _m(_cumprod)
2171
def _cumsum(self, axis=None, dtype=None, out=None):
2174
m = umath.logical_or.accumulate(self.mask, axis)
2175
return MaskedArray(data=self.filled(0).cumsum(axis, dtype), mask=m)
2176
array.cumsum = _m(_cumsum)
2178
array.diagonal = _m(diagonal)
2179
array.dump = _m(not_implemented)
2180
array.dumps = _m(not_implemented)
2181
array.fill = _m(not_implemented)
2182
array.flags = property(_m(not_implemented))
2183
array.flatten = _m(ravel)
2184
array.getfield = _m(not_implemented)
2186
def _max(a, axis=None, out=None):
2188
raise TypeError("Output arrays Unsupported for masked arrays")
2192
return maximum.reduce(a, axis)
2193
array.max = _m(_max)
2195
def _min(a, axis=None, out=None):
2197
raise TypeError("Output arrays Unsupported for masked arrays")
2201
return minimum.reduce(a, axis)
2202
array.min = _m(_min)
2204
array.mean = _m(new_average)
2205
array.nbytes = property(_m(not_implemented))
2206
array.newbyteorder = _m(not_implemented)
2207
array.nonzero = _m(nonzero)
2208
array.prod = _m(product)
2210
def _ptp(a,axis=None,out=None):
2211
return a.max(axis,out)-a.min(axis)
2212
array.ptp = _m(_ptp)
2213
array.repeat = _m(new_repeat)
2214
array.resize = _m(resize)
2215
array.searchsorted = _m(not_implemented)
2216
array.setfield = _m(not_implemented)
2217
array.setflags = _m(not_implemented)
2218
array.sort = _m(not_implemented) # NB: ndarray.sort is inplace
2222
result = MaskedArray(data = self.data.squeeze(),
2223
mask = self.mask.squeeze())
2224
except AttributeError:
2225
result = _wrapit(self, 'squeeze')
2227
array.squeeze = _m(_squeeze)
2229
array.strides = property(_m(not_implemented))
2231
def _swapaxes(self,axis1,axis2):
2232
return MaskedArray(data = self.data.swapaxes(axis1, axis2),
2233
mask = self.mask.swapaxes(axis1, axis2))
2234
array.swapaxes = _m(_swapaxes)
2235
array.take = _m(new_take)
2236
array.tofile = _m(not_implemented)
2237
array.trace = _m(trace)
2238
array.transpose = _m(transpose)
2240
def _var(self,axis=None,dtype=None, out=None):
2242
return numeric.asarray(self.compressed()).var()
2243
a = self.swapaxes(axis,0)
2244
a = a - a.mean(axis=0)
2246
a /= a.count(axis=0)
2247
return a.swapaxes(0,axis).sum(axis)
2248
def _std(self,axis=None, dtype=None, out=None):
2249
return (self.var(axis,dtype))**0.5
2250
array.var = _m(_var)
2251
array.std = _m(_std)
2253
array.view = _m(not_implemented)
2254
array.round = _m(around)
2255
del _m, MethodType, not_implemented
2258
masked = MaskedArray(0, int, mask=1)
7
2260
def repeat(a, repeats, axis=0):
8
return nca.repeat(a, repeats, axis)
2261
return new_repeat(a, repeats, axis)
10
2263
def average(a, axis=0, weights=None, returned=0):
11
return nca.average(a, axis, weights, returned)
2264
return new_average(a, axis, weights, returned)
13
2266
def take(a, indices, axis=0):
14
return nca.average(a, indices, axis=0)
2267
return new_take(a, indices, axis)