Chris@87
|
1 from __future__ import division, absolute_import, print_function
|
Chris@87
|
2
|
Chris@87
|
3 import sys
|
Chris@87
|
4 import math
|
Chris@87
|
5
|
Chris@87
|
6 import numpy.core.numeric as _nx
|
Chris@87
|
7 from numpy.core.numeric import (
|
Chris@87
|
8 asarray, ScalarType, array, alltrue, cumprod, arange
|
Chris@87
|
9 )
|
Chris@87
|
10 from numpy.core.numerictypes import find_common_type
|
Chris@87
|
11
|
Chris@87
|
12 from . import function_base
|
Chris@87
|
13 import numpy.matrixlib as matrix
|
Chris@87
|
14 from .function_base import diff
|
Chris@87
|
15 from numpy.lib._compiled_base import ravel_multi_index, unravel_index
|
Chris@87
|
16 from numpy.lib.stride_tricks import as_strided
|
Chris@87
|
17
|
Chris@87
|
18 makemat = matrix.matrix
|
Chris@87
|
19
|
Chris@87
|
20
|
Chris@87
|
21 __all__ = [
|
Chris@87
|
22 'ravel_multi_index', 'unravel_index', 'mgrid', 'ogrid', 'r_', 'c_',
|
Chris@87
|
23 's_', 'index_exp', 'ix_', 'ndenumerate', 'ndindex', 'fill_diagonal',
|
Chris@87
|
24 'diag_indices', 'diag_indices_from'
|
Chris@87
|
25 ]
|
Chris@87
|
26
|
Chris@87
|
27
|
Chris@87
|
28 def ix_(*args):
|
Chris@87
|
29 """
|
Chris@87
|
30 Construct an open mesh from multiple sequences.
|
Chris@87
|
31
|
Chris@87
|
32 This function takes N 1-D sequences and returns N outputs with N
|
Chris@87
|
33 dimensions each, such that the shape is 1 in all but one dimension
|
Chris@87
|
34 and the dimension with the non-unit shape value cycles through all
|
Chris@87
|
35 N dimensions.
|
Chris@87
|
36
|
Chris@87
|
37 Using `ix_` one can quickly construct index arrays that will index
|
Chris@87
|
38 the cross product. ``a[np.ix_([1,3],[2,5])]`` returns the array
|
Chris@87
|
39 ``[[a[1,2] a[1,5]], [a[3,2] a[3,5]]]``.
|
Chris@87
|
40
|
Chris@87
|
41 Parameters
|
Chris@87
|
42 ----------
|
Chris@87
|
43 args : 1-D sequences
|
Chris@87
|
44
|
Chris@87
|
45 Returns
|
Chris@87
|
46 -------
|
Chris@87
|
47 out : tuple of ndarrays
|
Chris@87
|
48 N arrays with N dimensions each, with N the number of input
|
Chris@87
|
49 sequences. Together these arrays form an open mesh.
|
Chris@87
|
50
|
Chris@87
|
51 See Also
|
Chris@87
|
52 --------
|
Chris@87
|
53 ogrid, mgrid, meshgrid
|
Chris@87
|
54
|
Chris@87
|
55 Examples
|
Chris@87
|
56 --------
|
Chris@87
|
57 >>> a = np.arange(10).reshape(2, 5)
|
Chris@87
|
58 >>> a
|
Chris@87
|
59 array([[0, 1, 2, 3, 4],
|
Chris@87
|
60 [5, 6, 7, 8, 9]])
|
Chris@87
|
61 >>> ixgrid = np.ix_([0,1], [2,4])
|
Chris@87
|
62 >>> ixgrid
|
Chris@87
|
63 (array([[0],
|
Chris@87
|
64 [1]]), array([[2, 4]]))
|
Chris@87
|
65 >>> ixgrid[0].shape, ixgrid[1].shape
|
Chris@87
|
66 ((2, 1), (1, 2))
|
Chris@87
|
67 >>> a[ixgrid]
|
Chris@87
|
68 array([[2, 4],
|
Chris@87
|
69 [7, 9]])
|
Chris@87
|
70
|
Chris@87
|
71 """
|
Chris@87
|
72 out = []
|
Chris@87
|
73 nd = len(args)
|
Chris@87
|
74 baseshape = [1]*nd
|
Chris@87
|
75 for k in range(nd):
|
Chris@87
|
76 new = _nx.asarray(args[k])
|
Chris@87
|
77 if (new.ndim != 1):
|
Chris@87
|
78 raise ValueError("Cross index must be 1 dimensional")
|
Chris@87
|
79 if issubclass(new.dtype.type, _nx.bool_):
|
Chris@87
|
80 new = new.nonzero()[0]
|
Chris@87
|
81 baseshape[k] = len(new)
|
Chris@87
|
82 new = new.reshape(tuple(baseshape))
|
Chris@87
|
83 out.append(new)
|
Chris@87
|
84 baseshape[k] = 1
|
Chris@87
|
85 return tuple(out)
|
Chris@87
|
86
|
Chris@87
|
87 class nd_grid(object):
|
Chris@87
|
88 """
|
Chris@87
|
89 Construct a multi-dimensional "meshgrid".
|
Chris@87
|
90
|
Chris@87
|
91 ``grid = nd_grid()`` creates an instance which will return a mesh-grid
|
Chris@87
|
92 when indexed. The dimension and number of the output arrays are equal
|
Chris@87
|
93 to the number of indexing dimensions. If the step length is not a
|
Chris@87
|
94 complex number, then the stop is not inclusive.
|
Chris@87
|
95
|
Chris@87
|
96 However, if the step length is a **complex number** (e.g. 5j), then the
|
Chris@87
|
97 integer part of its magnitude is interpreted as specifying the
|
Chris@87
|
98 number of points to create between the start and stop values, where
|
Chris@87
|
99 the stop value **is inclusive**.
|
Chris@87
|
100
|
Chris@87
|
101 If instantiated with an argument of ``sparse=True``, the mesh-grid is
|
Chris@87
|
102 open (or not fleshed out) so that only one-dimension of each returned
|
Chris@87
|
103 argument is greater than 1.
|
Chris@87
|
104
|
Chris@87
|
105 Parameters
|
Chris@87
|
106 ----------
|
Chris@87
|
107 sparse : bool, optional
|
Chris@87
|
108 Whether the grid is sparse or not. Default is False.
|
Chris@87
|
109
|
Chris@87
|
110 Notes
|
Chris@87
|
111 -----
|
Chris@87
|
112 Two instances of `nd_grid` are made available in the NumPy namespace,
|
Chris@87
|
113 `mgrid` and `ogrid`::
|
Chris@87
|
114
|
Chris@87
|
115 mgrid = nd_grid(sparse=False)
|
Chris@87
|
116 ogrid = nd_grid(sparse=True)
|
Chris@87
|
117
|
Chris@87
|
118 Users should use these pre-defined instances instead of using `nd_grid`
|
Chris@87
|
119 directly.
|
Chris@87
|
120
|
Chris@87
|
121 Examples
|
Chris@87
|
122 --------
|
Chris@87
|
123 >>> mgrid = np.lib.index_tricks.nd_grid()
|
Chris@87
|
124 >>> mgrid[0:5,0:5]
|
Chris@87
|
125 array([[[0, 0, 0, 0, 0],
|
Chris@87
|
126 [1, 1, 1, 1, 1],
|
Chris@87
|
127 [2, 2, 2, 2, 2],
|
Chris@87
|
128 [3, 3, 3, 3, 3],
|
Chris@87
|
129 [4, 4, 4, 4, 4]],
|
Chris@87
|
130 [[0, 1, 2, 3, 4],
|
Chris@87
|
131 [0, 1, 2, 3, 4],
|
Chris@87
|
132 [0, 1, 2, 3, 4],
|
Chris@87
|
133 [0, 1, 2, 3, 4],
|
Chris@87
|
134 [0, 1, 2, 3, 4]]])
|
Chris@87
|
135 >>> mgrid[-1:1:5j]
|
Chris@87
|
136 array([-1. , -0.5, 0. , 0.5, 1. ])
|
Chris@87
|
137
|
Chris@87
|
138 >>> ogrid = np.lib.index_tricks.nd_grid(sparse=True)
|
Chris@87
|
139 >>> ogrid[0:5,0:5]
|
Chris@87
|
140 [array([[0],
|
Chris@87
|
141 [1],
|
Chris@87
|
142 [2],
|
Chris@87
|
143 [3],
|
Chris@87
|
144 [4]]), array([[0, 1, 2, 3, 4]])]
|
Chris@87
|
145
|
Chris@87
|
146 """
|
Chris@87
|
147
|
Chris@87
|
148 def __init__(self, sparse=False):
|
Chris@87
|
149 self.sparse = sparse
|
Chris@87
|
150
|
Chris@87
|
151 def __getitem__(self, key):
|
Chris@87
|
152 try:
|
Chris@87
|
153 size = []
|
Chris@87
|
154 typ = int
|
Chris@87
|
155 for k in range(len(key)):
|
Chris@87
|
156 step = key[k].step
|
Chris@87
|
157 start = key[k].start
|
Chris@87
|
158 if start is None:
|
Chris@87
|
159 start = 0
|
Chris@87
|
160 if step is None:
|
Chris@87
|
161 step = 1
|
Chris@87
|
162 if isinstance(step, complex):
|
Chris@87
|
163 size.append(int(abs(step)))
|
Chris@87
|
164 typ = float
|
Chris@87
|
165 else:
|
Chris@87
|
166 size.append(
|
Chris@87
|
167 int(math.ceil((key[k].stop - start)/(step*1.0))))
|
Chris@87
|
168 if (isinstance(step, float) or
|
Chris@87
|
169 isinstance(start, float) or
|
Chris@87
|
170 isinstance(key[k].stop, float)):
|
Chris@87
|
171 typ = float
|
Chris@87
|
172 if self.sparse:
|
Chris@87
|
173 nn = [_nx.arange(_x, dtype=_t)
|
Chris@87
|
174 for _x, _t in zip(size, (typ,)*len(size))]
|
Chris@87
|
175 else:
|
Chris@87
|
176 nn = _nx.indices(size, typ)
|
Chris@87
|
177 for k in range(len(size)):
|
Chris@87
|
178 step = key[k].step
|
Chris@87
|
179 start = key[k].start
|
Chris@87
|
180 if start is None:
|
Chris@87
|
181 start = 0
|
Chris@87
|
182 if step is None:
|
Chris@87
|
183 step = 1
|
Chris@87
|
184 if isinstance(step, complex):
|
Chris@87
|
185 step = int(abs(step))
|
Chris@87
|
186 if step != 1:
|
Chris@87
|
187 step = (key[k].stop - start)/float(step-1)
|
Chris@87
|
188 nn[k] = (nn[k]*step+start)
|
Chris@87
|
189 if self.sparse:
|
Chris@87
|
190 slobj = [_nx.newaxis]*len(size)
|
Chris@87
|
191 for k in range(len(size)):
|
Chris@87
|
192 slobj[k] = slice(None, None)
|
Chris@87
|
193 nn[k] = nn[k][slobj]
|
Chris@87
|
194 slobj[k] = _nx.newaxis
|
Chris@87
|
195 return nn
|
Chris@87
|
196 except (IndexError, TypeError):
|
Chris@87
|
197 step = key.step
|
Chris@87
|
198 stop = key.stop
|
Chris@87
|
199 start = key.start
|
Chris@87
|
200 if start is None:
|
Chris@87
|
201 start = 0
|
Chris@87
|
202 if isinstance(step, complex):
|
Chris@87
|
203 step = abs(step)
|
Chris@87
|
204 length = int(step)
|
Chris@87
|
205 if step != 1:
|
Chris@87
|
206 step = (key.stop-start)/float(step-1)
|
Chris@87
|
207 stop = key.stop + step
|
Chris@87
|
208 return _nx.arange(0, length, 1, float)*step + start
|
Chris@87
|
209 else:
|
Chris@87
|
210 return _nx.arange(start, stop, step)
|
Chris@87
|
211
|
Chris@87
|
212 def __getslice__(self, i, j):
|
Chris@87
|
213 return _nx.arange(i, j)
|
Chris@87
|
214
|
Chris@87
|
215 def __len__(self):
|
Chris@87
|
216 return 0
|
Chris@87
|
217
|
Chris@87
|
218 mgrid = nd_grid(sparse=False)
|
Chris@87
|
219 ogrid = nd_grid(sparse=True)
|
Chris@87
|
220 mgrid.__doc__ = None # set in numpy.add_newdocs
|
Chris@87
|
221 ogrid.__doc__ = None # set in numpy.add_newdocs
|
Chris@87
|
222
|
Chris@87
|
223 class AxisConcatenator(object):
|
Chris@87
|
224 """
|
Chris@87
|
225 Translates slice objects to concatenation along an axis.
|
Chris@87
|
226
|
Chris@87
|
227 For detailed documentation on usage, see `r_`.
|
Chris@87
|
228
|
Chris@87
|
229 """
|
Chris@87
|
230
|
Chris@87
|
231 def _retval(self, res):
|
Chris@87
|
232 if self.matrix:
|
Chris@87
|
233 oldndim = res.ndim
|
Chris@87
|
234 res = makemat(res)
|
Chris@87
|
235 if oldndim == 1 and self.col:
|
Chris@87
|
236 res = res.T
|
Chris@87
|
237 self.axis = self._axis
|
Chris@87
|
238 self.matrix = self._matrix
|
Chris@87
|
239 self.col = 0
|
Chris@87
|
240 return res
|
Chris@87
|
241
|
Chris@87
|
242 def __init__(self, axis=0, matrix=False, ndmin=1, trans1d=-1):
|
Chris@87
|
243 self._axis = axis
|
Chris@87
|
244 self._matrix = matrix
|
Chris@87
|
245 self.axis = axis
|
Chris@87
|
246 self.matrix = matrix
|
Chris@87
|
247 self.col = 0
|
Chris@87
|
248 self.trans1d = trans1d
|
Chris@87
|
249 self.ndmin = ndmin
|
Chris@87
|
250
|
Chris@87
|
251 def __getitem__(self, key):
|
Chris@87
|
252 trans1d = self.trans1d
|
Chris@87
|
253 ndmin = self.ndmin
|
Chris@87
|
254 if isinstance(key, str):
|
Chris@87
|
255 frame = sys._getframe().f_back
|
Chris@87
|
256 mymat = matrix.bmat(key, frame.f_globals, frame.f_locals)
|
Chris@87
|
257 return mymat
|
Chris@87
|
258 if not isinstance(key, tuple):
|
Chris@87
|
259 key = (key,)
|
Chris@87
|
260 objs = []
|
Chris@87
|
261 scalars = []
|
Chris@87
|
262 arraytypes = []
|
Chris@87
|
263 scalartypes = []
|
Chris@87
|
264 for k in range(len(key)):
|
Chris@87
|
265 scalar = False
|
Chris@87
|
266 if isinstance(key[k], slice):
|
Chris@87
|
267 step = key[k].step
|
Chris@87
|
268 start = key[k].start
|
Chris@87
|
269 stop = key[k].stop
|
Chris@87
|
270 if start is None:
|
Chris@87
|
271 start = 0
|
Chris@87
|
272 if step is None:
|
Chris@87
|
273 step = 1
|
Chris@87
|
274 if isinstance(step, complex):
|
Chris@87
|
275 size = int(abs(step))
|
Chris@87
|
276 newobj = function_base.linspace(start, stop, num=size)
|
Chris@87
|
277 else:
|
Chris@87
|
278 newobj = _nx.arange(start, stop, step)
|
Chris@87
|
279 if ndmin > 1:
|
Chris@87
|
280 newobj = array(newobj, copy=False, ndmin=ndmin)
|
Chris@87
|
281 if trans1d != -1:
|
Chris@87
|
282 newobj = newobj.swapaxes(-1, trans1d)
|
Chris@87
|
283 elif isinstance(key[k], str):
|
Chris@87
|
284 if k != 0:
|
Chris@87
|
285 raise ValueError("special directives must be the "
|
Chris@87
|
286 "first entry.")
|
Chris@87
|
287 key0 = key[0]
|
Chris@87
|
288 if key0 in 'rc':
|
Chris@87
|
289 self.matrix = True
|
Chris@87
|
290 self.col = (key0 == 'c')
|
Chris@87
|
291 continue
|
Chris@87
|
292 if ',' in key0:
|
Chris@87
|
293 vec = key0.split(',')
|
Chris@87
|
294 try:
|
Chris@87
|
295 self.axis, ndmin = \
|
Chris@87
|
296 [int(x) for x in vec[:2]]
|
Chris@87
|
297 if len(vec) == 3:
|
Chris@87
|
298 trans1d = int(vec[2])
|
Chris@87
|
299 continue
|
Chris@87
|
300 except:
|
Chris@87
|
301 raise ValueError("unknown special directive")
|
Chris@87
|
302 try:
|
Chris@87
|
303 self.axis = int(key[k])
|
Chris@87
|
304 continue
|
Chris@87
|
305 except (ValueError, TypeError):
|
Chris@87
|
306 raise ValueError("unknown special directive")
|
Chris@87
|
307 elif type(key[k]) in ScalarType:
|
Chris@87
|
308 newobj = array(key[k], ndmin=ndmin)
|
Chris@87
|
309 scalars.append(k)
|
Chris@87
|
310 scalar = True
|
Chris@87
|
311 scalartypes.append(newobj.dtype)
|
Chris@87
|
312 else:
|
Chris@87
|
313 newobj = key[k]
|
Chris@87
|
314 if ndmin > 1:
|
Chris@87
|
315 tempobj = array(newobj, copy=False, subok=True)
|
Chris@87
|
316 newobj = array(newobj, copy=False, subok=True,
|
Chris@87
|
317 ndmin=ndmin)
|
Chris@87
|
318 if trans1d != -1 and tempobj.ndim < ndmin:
|
Chris@87
|
319 k2 = ndmin-tempobj.ndim
|
Chris@87
|
320 if (trans1d < 0):
|
Chris@87
|
321 trans1d += k2 + 1
|
Chris@87
|
322 defaxes = list(range(ndmin))
|
Chris@87
|
323 k1 = trans1d
|
Chris@87
|
324 axes = defaxes[:k1] + defaxes[k2:] + \
|
Chris@87
|
325 defaxes[k1:k2]
|
Chris@87
|
326 newobj = newobj.transpose(axes)
|
Chris@87
|
327 del tempobj
|
Chris@87
|
328 objs.append(newobj)
|
Chris@87
|
329 if not scalar and isinstance(newobj, _nx.ndarray):
|
Chris@87
|
330 arraytypes.append(newobj.dtype)
|
Chris@87
|
331
|
Chris@87
|
332 # Esure that scalars won't up-cast unless warranted
|
Chris@87
|
333 final_dtype = find_common_type(arraytypes, scalartypes)
|
Chris@87
|
334 if final_dtype is not None:
|
Chris@87
|
335 for k in scalars:
|
Chris@87
|
336 objs[k] = objs[k].astype(final_dtype)
|
Chris@87
|
337
|
Chris@87
|
338 res = _nx.concatenate(tuple(objs), axis=self.axis)
|
Chris@87
|
339 return self._retval(res)
|
Chris@87
|
340
|
Chris@87
|
341 def __getslice__(self, i, j):
|
Chris@87
|
342 res = _nx.arange(i, j)
|
Chris@87
|
343 return self._retval(res)
|
Chris@87
|
344
|
Chris@87
|
345 def __len__(self):
|
Chris@87
|
346 return 0
|
Chris@87
|
347
|
Chris@87
|
348 # separate classes are used here instead of just making r_ = concatentor(0),
|
Chris@87
|
349 # etc. because otherwise we couldn't get the doc string to come out right
|
Chris@87
|
350 # in help(r_)
|
Chris@87
|
351
|
Chris@87
|
352 class RClass(AxisConcatenator):
|
Chris@87
|
353 """
|
Chris@87
|
354 Translates slice objects to concatenation along the first axis.
|
Chris@87
|
355
|
Chris@87
|
356 This is a simple way to build up arrays quickly. There are two use cases.
|
Chris@87
|
357
|
Chris@87
|
358 1. If the index expression contains comma separated arrays, then stack
|
Chris@87
|
359 them along their first axis.
|
Chris@87
|
360 2. If the index expression contains slice notation or scalars then create
|
Chris@87
|
361 a 1-D array with a range indicated by the slice notation.
|
Chris@87
|
362
|
Chris@87
|
363 If slice notation is used, the syntax ``start:stop:step`` is equivalent
|
Chris@87
|
364 to ``np.arange(start, stop, step)`` inside of the brackets. However, if
|
Chris@87
|
365 ``step`` is an imaginary number (i.e. 100j) then its integer portion is
|
Chris@87
|
366 interpreted as a number-of-points desired and the start and stop are
|
Chris@87
|
367 inclusive. In other words ``start:stop:stepj`` is interpreted as
|
Chris@87
|
368 ``np.linspace(start, stop, step, endpoint=1)`` inside of the brackets.
|
Chris@87
|
369 After expansion of slice notation, all comma separated sequences are
|
Chris@87
|
370 concatenated together.
|
Chris@87
|
371
|
Chris@87
|
372 Optional character strings placed as the first element of the index
|
Chris@87
|
373 expression can be used to change the output. The strings 'r' or 'c' result
|
Chris@87
|
374 in matrix output. If the result is 1-D and 'r' is specified a 1 x N (row)
|
Chris@87
|
375 matrix is produced. If the result is 1-D and 'c' is specified, then a N x 1
|
Chris@87
|
376 (column) matrix is produced. If the result is 2-D then both provide the
|
Chris@87
|
377 same matrix result.
|
Chris@87
|
378
|
Chris@87
|
379 A string integer specifies which axis to stack multiple comma separated
|
Chris@87
|
380 arrays along. A string of two comma-separated integers allows indication
|
Chris@87
|
381 of the minimum number of dimensions to force each entry into as the
|
Chris@87
|
382 second integer (the axis to concatenate along is still the first integer).
|
Chris@87
|
383
|
Chris@87
|
384 A string with three comma-separated integers allows specification of the
|
Chris@87
|
385 axis to concatenate along, the minimum number of dimensions to force the
|
Chris@87
|
386 entries to, and which axis should contain the start of the arrays which
|
Chris@87
|
387 are less than the specified number of dimensions. In other words the third
|
Chris@87
|
388 integer allows you to specify where the 1's should be placed in the shape
|
Chris@87
|
389 of the arrays that have their shapes upgraded. By default, they are placed
|
Chris@87
|
390 in the front of the shape tuple. The third argument allows you to specify
|
Chris@87
|
391 where the start of the array should be instead. Thus, a third argument of
|
Chris@87
|
392 '0' would place the 1's at the end of the array shape. Negative integers
|
Chris@87
|
393 specify where in the new shape tuple the last dimension of upgraded arrays
|
Chris@87
|
394 should be placed, so the default is '-1'.
|
Chris@87
|
395
|
Chris@87
|
396 Parameters
|
Chris@87
|
397 ----------
|
Chris@87
|
398 Not a function, so takes no parameters
|
Chris@87
|
399
|
Chris@87
|
400
|
Chris@87
|
401 Returns
|
Chris@87
|
402 -------
|
Chris@87
|
403 A concatenated ndarray or matrix.
|
Chris@87
|
404
|
Chris@87
|
405 See Also
|
Chris@87
|
406 --------
|
Chris@87
|
407 concatenate : Join a sequence of arrays together.
|
Chris@87
|
408 c_ : Translates slice objects to concatenation along the second axis.
|
Chris@87
|
409
|
Chris@87
|
410 Examples
|
Chris@87
|
411 --------
|
Chris@87
|
412 >>> np.r_[np.array([1,2,3]), 0, 0, np.array([4,5,6])]
|
Chris@87
|
413 array([1, 2, 3, 0, 0, 4, 5, 6])
|
Chris@87
|
414 >>> np.r_[-1:1:6j, [0]*3, 5, 6]
|
Chris@87
|
415 array([-1. , -0.6, -0.2, 0.2, 0.6, 1. , 0. , 0. , 0. , 5. , 6. ])
|
Chris@87
|
416
|
Chris@87
|
417 String integers specify the axis to concatenate along or the minimum
|
Chris@87
|
418 number of dimensions to force entries into.
|
Chris@87
|
419
|
Chris@87
|
420 >>> a = np.array([[0, 1, 2], [3, 4, 5]])
|
Chris@87
|
421 >>> np.r_['-1', a, a] # concatenate along last axis
|
Chris@87
|
422 array([[0, 1, 2, 0, 1, 2],
|
Chris@87
|
423 [3, 4, 5, 3, 4, 5]])
|
Chris@87
|
424 >>> np.r_['0,2', [1,2,3], [4,5,6]] # concatenate along first axis, dim>=2
|
Chris@87
|
425 array([[1, 2, 3],
|
Chris@87
|
426 [4, 5, 6]])
|
Chris@87
|
427
|
Chris@87
|
428 >>> np.r_['0,2,0', [1,2,3], [4,5,6]]
|
Chris@87
|
429 array([[1],
|
Chris@87
|
430 [2],
|
Chris@87
|
431 [3],
|
Chris@87
|
432 [4],
|
Chris@87
|
433 [5],
|
Chris@87
|
434 [6]])
|
Chris@87
|
435 >>> np.r_['1,2,0', [1,2,3], [4,5,6]]
|
Chris@87
|
436 array([[1, 4],
|
Chris@87
|
437 [2, 5],
|
Chris@87
|
438 [3, 6]])
|
Chris@87
|
439
|
Chris@87
|
440 Using 'r' or 'c' as a first string argument creates a matrix.
|
Chris@87
|
441
|
Chris@87
|
442 >>> np.r_['r',[1,2,3], [4,5,6]]
|
Chris@87
|
443 matrix([[1, 2, 3, 4, 5, 6]])
|
Chris@87
|
444
|
Chris@87
|
445 """
|
Chris@87
|
446
|
Chris@87
|
447 def __init__(self):
|
Chris@87
|
448 AxisConcatenator.__init__(self, 0)
|
Chris@87
|
449
|
Chris@87
|
450 r_ = RClass()
|
Chris@87
|
451
|
Chris@87
|
452 class CClass(AxisConcatenator):
|
Chris@87
|
453 """
|
Chris@87
|
454 Translates slice objects to concatenation along the second axis.
|
Chris@87
|
455
|
Chris@87
|
456 This is short-hand for ``np.r_['-1,2,0', index expression]``, which is
|
Chris@87
|
457 useful because of its common occurrence. In particular, arrays will be
|
Chris@87
|
458 stacked along their last axis after being upgraded to at least 2-D with
|
Chris@87
|
459 1's post-pended to the shape (column vectors made out of 1-D arrays).
|
Chris@87
|
460
|
Chris@87
|
461 For detailed documentation, see `r_`.
|
Chris@87
|
462
|
Chris@87
|
463 Examples
|
Chris@87
|
464 --------
|
Chris@87
|
465 >>> np.c_[np.array([[1,2,3]]), 0, 0, np.array([[4,5,6]])]
|
Chris@87
|
466 array([[1, 2, 3, 0, 0, 4, 5, 6]])
|
Chris@87
|
467
|
Chris@87
|
468 """
|
Chris@87
|
469
|
Chris@87
|
470 def __init__(self):
|
Chris@87
|
471 AxisConcatenator.__init__(self, -1, ndmin=2, trans1d=0)
|
Chris@87
|
472
|
Chris@87
|
473 c_ = CClass()
|
Chris@87
|
474
|
Chris@87
|
475 class ndenumerate(object):
|
Chris@87
|
476 """
|
Chris@87
|
477 Multidimensional index iterator.
|
Chris@87
|
478
|
Chris@87
|
479 Return an iterator yielding pairs of array coordinates and values.
|
Chris@87
|
480
|
Chris@87
|
481 Parameters
|
Chris@87
|
482 ----------
|
Chris@87
|
483 a : ndarray
|
Chris@87
|
484 Input array.
|
Chris@87
|
485
|
Chris@87
|
486 See Also
|
Chris@87
|
487 --------
|
Chris@87
|
488 ndindex, flatiter
|
Chris@87
|
489
|
Chris@87
|
490 Examples
|
Chris@87
|
491 --------
|
Chris@87
|
492 >>> a = np.array([[1, 2], [3, 4]])
|
Chris@87
|
493 >>> for index, x in np.ndenumerate(a):
|
Chris@87
|
494 ... print index, x
|
Chris@87
|
495 (0, 0) 1
|
Chris@87
|
496 (0, 1) 2
|
Chris@87
|
497 (1, 0) 3
|
Chris@87
|
498 (1, 1) 4
|
Chris@87
|
499
|
Chris@87
|
500 """
|
Chris@87
|
501
|
Chris@87
|
502 def __init__(self, arr):
|
Chris@87
|
503 self.iter = asarray(arr).flat
|
Chris@87
|
504
|
Chris@87
|
505 def __next__(self):
|
Chris@87
|
506 """
|
Chris@87
|
507 Standard iterator method, returns the index tuple and array value.
|
Chris@87
|
508
|
Chris@87
|
509 Returns
|
Chris@87
|
510 -------
|
Chris@87
|
511 coords : tuple of ints
|
Chris@87
|
512 The indices of the current iteration.
|
Chris@87
|
513 val : scalar
|
Chris@87
|
514 The array element of the current iteration.
|
Chris@87
|
515
|
Chris@87
|
516 """
|
Chris@87
|
517 return self.iter.coords, next(self.iter)
|
Chris@87
|
518
|
Chris@87
|
519 def __iter__(self):
|
Chris@87
|
520 return self
|
Chris@87
|
521
|
Chris@87
|
522 next = __next__
|
Chris@87
|
523
|
Chris@87
|
524
|
Chris@87
|
525 class ndindex(object):
|
Chris@87
|
526 """
|
Chris@87
|
527 An N-dimensional iterator object to index arrays.
|
Chris@87
|
528
|
Chris@87
|
529 Given the shape of an array, an `ndindex` instance iterates over
|
Chris@87
|
530 the N-dimensional index of the array. At each iteration a tuple
|
Chris@87
|
531 of indices is returned, the last dimension is iterated over first.
|
Chris@87
|
532
|
Chris@87
|
533 Parameters
|
Chris@87
|
534 ----------
|
Chris@87
|
535 `*args` : ints
|
Chris@87
|
536 The size of each dimension of the array.
|
Chris@87
|
537
|
Chris@87
|
538 See Also
|
Chris@87
|
539 --------
|
Chris@87
|
540 ndenumerate, flatiter
|
Chris@87
|
541
|
Chris@87
|
542 Examples
|
Chris@87
|
543 --------
|
Chris@87
|
544 >>> for index in np.ndindex(3, 2, 1):
|
Chris@87
|
545 ... print index
|
Chris@87
|
546 (0, 0, 0)
|
Chris@87
|
547 (0, 1, 0)
|
Chris@87
|
548 (1, 0, 0)
|
Chris@87
|
549 (1, 1, 0)
|
Chris@87
|
550 (2, 0, 0)
|
Chris@87
|
551 (2, 1, 0)
|
Chris@87
|
552
|
Chris@87
|
553 """
|
Chris@87
|
554
|
Chris@87
|
555 def __init__(self, *shape):
|
Chris@87
|
556 if len(shape) == 1 and isinstance(shape[0], tuple):
|
Chris@87
|
557 shape = shape[0]
|
Chris@87
|
558 x = as_strided(_nx.zeros(1), shape=shape,
|
Chris@87
|
559 strides=_nx.zeros_like(shape))
|
Chris@87
|
560 self._it = _nx.nditer(x, flags=['multi_index', 'zerosize_ok'],
|
Chris@87
|
561 order='C')
|
Chris@87
|
562
|
Chris@87
|
563 def __iter__(self):
|
Chris@87
|
564 return self
|
Chris@87
|
565
|
Chris@87
|
566 def ndincr(self):
|
Chris@87
|
567 """
|
Chris@87
|
568 Increment the multi-dimensional index by one.
|
Chris@87
|
569
|
Chris@87
|
570 This method is for backward compatibility only: do not use.
|
Chris@87
|
571 """
|
Chris@87
|
572 next(self)
|
Chris@87
|
573
|
Chris@87
|
574 def __next__(self):
|
Chris@87
|
575 """
|
Chris@87
|
576 Standard iterator method, updates the index and returns the index
|
Chris@87
|
577 tuple.
|
Chris@87
|
578
|
Chris@87
|
579 Returns
|
Chris@87
|
580 -------
|
Chris@87
|
581 val : tuple of ints
|
Chris@87
|
582 Returns a tuple containing the indices of the current
|
Chris@87
|
583 iteration.
|
Chris@87
|
584
|
Chris@87
|
585 """
|
Chris@87
|
586 next(self._it)
|
Chris@87
|
587 return self._it.multi_index
|
Chris@87
|
588
|
Chris@87
|
589 next = __next__
|
Chris@87
|
590
|
Chris@87
|
591
|
Chris@87
|
592 # You can do all this with slice() plus a few special objects,
|
Chris@87
|
593 # but there's a lot to remember. This version is simpler because
|
Chris@87
|
594 # it uses the standard array indexing syntax.
|
Chris@87
|
595 #
|
Chris@87
|
596 # Written by Konrad Hinsen <hinsen@cnrs-orleans.fr>
|
Chris@87
|
597 # last revision: 1999-7-23
|
Chris@87
|
598 #
|
Chris@87
|
599 # Cosmetic changes by T. Oliphant 2001
|
Chris@87
|
600 #
|
Chris@87
|
601 #
|
Chris@87
|
602
|
Chris@87
|
603 class IndexExpression(object):
|
Chris@87
|
604 """
|
Chris@87
|
605 A nicer way to build up index tuples for arrays.
|
Chris@87
|
606
|
Chris@87
|
607 .. note::
|
Chris@87
|
608 Use one of the two predefined instances `index_exp` or `s_`
|
Chris@87
|
609 rather than directly using `IndexExpression`.
|
Chris@87
|
610
|
Chris@87
|
611 For any index combination, including slicing and axis insertion,
|
Chris@87
|
612 ``a[indices]`` is the same as ``a[np.index_exp[indices]]`` for any
|
Chris@87
|
613 array `a`. However, ``np.index_exp[indices]`` can be used anywhere
|
Chris@87
|
614 in Python code and returns a tuple of slice objects that can be
|
Chris@87
|
615 used in the construction of complex index expressions.
|
Chris@87
|
616
|
Chris@87
|
617 Parameters
|
Chris@87
|
618 ----------
|
Chris@87
|
619 maketuple : bool
|
Chris@87
|
620 If True, always returns a tuple.
|
Chris@87
|
621
|
Chris@87
|
622 See Also
|
Chris@87
|
623 --------
|
Chris@87
|
624 index_exp : Predefined instance that always returns a tuple:
|
Chris@87
|
625 `index_exp = IndexExpression(maketuple=True)`.
|
Chris@87
|
626 s_ : Predefined instance without tuple conversion:
|
Chris@87
|
627 `s_ = IndexExpression(maketuple=False)`.
|
Chris@87
|
628
|
Chris@87
|
629 Notes
|
Chris@87
|
630 -----
|
Chris@87
|
631 You can do all this with `slice()` plus a few special objects,
|
Chris@87
|
632 but there's a lot to remember and this version is simpler because
|
Chris@87
|
633 it uses the standard array indexing syntax.
|
Chris@87
|
634
|
Chris@87
|
635 Examples
|
Chris@87
|
636 --------
|
Chris@87
|
637 >>> np.s_[2::2]
|
Chris@87
|
638 slice(2, None, 2)
|
Chris@87
|
639 >>> np.index_exp[2::2]
|
Chris@87
|
640 (slice(2, None, 2),)
|
Chris@87
|
641
|
Chris@87
|
642 >>> np.array([0, 1, 2, 3, 4])[np.s_[2::2]]
|
Chris@87
|
643 array([2, 4])
|
Chris@87
|
644
|
Chris@87
|
645 """
|
Chris@87
|
646
|
Chris@87
|
647 def __init__(self, maketuple):
|
Chris@87
|
648 self.maketuple = maketuple
|
Chris@87
|
649
|
Chris@87
|
650 def __getitem__(self, item):
|
Chris@87
|
651 if self.maketuple and not isinstance(item, tuple):
|
Chris@87
|
652 return (item,)
|
Chris@87
|
653 else:
|
Chris@87
|
654 return item
|
Chris@87
|
655
|
Chris@87
|
656 index_exp = IndexExpression(maketuple=True)
|
Chris@87
|
657 s_ = IndexExpression(maketuple=False)
|
Chris@87
|
658
|
Chris@87
|
659 # End contribution from Konrad.
|
Chris@87
|
660
|
Chris@87
|
661
|
Chris@87
|
662 # The following functions complement those in twodim_base, but are
|
Chris@87
|
663 # applicable to N-dimensions.
|
Chris@87
|
664
|
Chris@87
|
665 def fill_diagonal(a, val, wrap=False):
|
Chris@87
|
666 """Fill the main diagonal of the given array of any dimensionality.
|
Chris@87
|
667
|
Chris@87
|
668 For an array `a` with ``a.ndim > 2``, the diagonal is the list of
|
Chris@87
|
669 locations with indices ``a[i, i, ..., i]`` all identical. This function
|
Chris@87
|
670 modifies the input array in-place, it does not return a value.
|
Chris@87
|
671
|
Chris@87
|
672 Parameters
|
Chris@87
|
673 ----------
|
Chris@87
|
674 a : array, at least 2-D.
|
Chris@87
|
675 Array whose diagonal is to be filled, it gets modified in-place.
|
Chris@87
|
676
|
Chris@87
|
677 val : scalar
|
Chris@87
|
678 Value to be written on the diagonal, its type must be compatible with
|
Chris@87
|
679 that of the array a.
|
Chris@87
|
680
|
Chris@87
|
681 wrap : bool
|
Chris@87
|
682 For tall matrices in NumPy version up to 1.6.2, the
|
Chris@87
|
683 diagonal "wrapped" after N columns. You can have this behavior
|
Chris@87
|
684 with this option. This affect only tall matrices.
|
Chris@87
|
685
|
Chris@87
|
686 See also
|
Chris@87
|
687 --------
|
Chris@87
|
688 diag_indices, diag_indices_from
|
Chris@87
|
689
|
Chris@87
|
690 Notes
|
Chris@87
|
691 -----
|
Chris@87
|
692 .. versionadded:: 1.4.0
|
Chris@87
|
693
|
Chris@87
|
694 This functionality can be obtained via `diag_indices`, but internally
|
Chris@87
|
695 this version uses a much faster implementation that never constructs the
|
Chris@87
|
696 indices and uses simple slicing.
|
Chris@87
|
697
|
Chris@87
|
698 Examples
|
Chris@87
|
699 --------
|
Chris@87
|
700 >>> a = np.zeros((3, 3), int)
|
Chris@87
|
701 >>> np.fill_diagonal(a, 5)
|
Chris@87
|
702 >>> a
|
Chris@87
|
703 array([[5, 0, 0],
|
Chris@87
|
704 [0, 5, 0],
|
Chris@87
|
705 [0, 0, 5]])
|
Chris@87
|
706
|
Chris@87
|
707 The same function can operate on a 4-D array:
|
Chris@87
|
708
|
Chris@87
|
709 >>> a = np.zeros((3, 3, 3, 3), int)
|
Chris@87
|
710 >>> np.fill_diagonal(a, 4)
|
Chris@87
|
711
|
Chris@87
|
712 We only show a few blocks for clarity:
|
Chris@87
|
713
|
Chris@87
|
714 >>> a[0, 0]
|
Chris@87
|
715 array([[4, 0, 0],
|
Chris@87
|
716 [0, 0, 0],
|
Chris@87
|
717 [0, 0, 0]])
|
Chris@87
|
718 >>> a[1, 1]
|
Chris@87
|
719 array([[0, 0, 0],
|
Chris@87
|
720 [0, 4, 0],
|
Chris@87
|
721 [0, 0, 0]])
|
Chris@87
|
722 >>> a[2, 2]
|
Chris@87
|
723 array([[0, 0, 0],
|
Chris@87
|
724 [0, 0, 0],
|
Chris@87
|
725 [0, 0, 4]])
|
Chris@87
|
726
|
Chris@87
|
727 # tall matrices no wrap
|
Chris@87
|
728 >>> a = np.zeros((5, 3),int)
|
Chris@87
|
729 >>> fill_diagonal(a, 4)
|
Chris@87
|
730 array([[4, 0, 0],
|
Chris@87
|
731 [0, 4, 0],
|
Chris@87
|
732 [0, 0, 4],
|
Chris@87
|
733 [0, 0, 0],
|
Chris@87
|
734 [0, 0, 0]])
|
Chris@87
|
735
|
Chris@87
|
736 # tall matrices wrap
|
Chris@87
|
737 >>> a = np.zeros((5, 3),int)
|
Chris@87
|
738 >>> fill_diagonal(a, 4)
|
Chris@87
|
739 array([[4, 0, 0],
|
Chris@87
|
740 [0, 4, 0],
|
Chris@87
|
741 [0, 0, 4],
|
Chris@87
|
742 [0, 0, 0],
|
Chris@87
|
743 [4, 0, 0]])
|
Chris@87
|
744
|
Chris@87
|
745 # wide matrices
|
Chris@87
|
746 >>> a = np.zeros((3, 5),int)
|
Chris@87
|
747 >>> fill_diagonal(a, 4)
|
Chris@87
|
748 array([[4, 0, 0, 0, 0],
|
Chris@87
|
749 [0, 4, 0, 0, 0],
|
Chris@87
|
750 [0, 0, 4, 0, 0]])
|
Chris@87
|
751
|
Chris@87
|
752 """
|
Chris@87
|
753 if a.ndim < 2:
|
Chris@87
|
754 raise ValueError("array must be at least 2-d")
|
Chris@87
|
755 end = None
|
Chris@87
|
756 if a.ndim == 2:
|
Chris@87
|
757 # Explicit, fast formula for the common case. For 2-d arrays, we
|
Chris@87
|
758 # accept rectangular ones.
|
Chris@87
|
759 step = a.shape[1] + 1
|
Chris@87
|
760 #This is needed to don't have tall matrix have the diagonal wrap.
|
Chris@87
|
761 if not wrap:
|
Chris@87
|
762 end = a.shape[1] * a.shape[1]
|
Chris@87
|
763 else:
|
Chris@87
|
764 # For more than d=2, the strided formula is only valid for arrays with
|
Chris@87
|
765 # all dimensions equal, so we check first.
|
Chris@87
|
766 if not alltrue(diff(a.shape) == 0):
|
Chris@87
|
767 raise ValueError("All dimensions of input must be of equal length")
|
Chris@87
|
768 step = 1 + (cumprod(a.shape[:-1])).sum()
|
Chris@87
|
769
|
Chris@87
|
770 # Write the value out into the diagonal.
|
Chris@87
|
771 a.flat[:end:step] = val
|
Chris@87
|
772
|
Chris@87
|
773
|
Chris@87
|
774 def diag_indices(n, ndim=2):
|
Chris@87
|
775 """
|
Chris@87
|
776 Return the indices to access the main diagonal of an array.
|
Chris@87
|
777
|
Chris@87
|
778 This returns a tuple of indices that can be used to access the main
|
Chris@87
|
779 diagonal of an array `a` with ``a.ndim >= 2`` dimensions and shape
|
Chris@87
|
780 (n, n, ..., n). For ``a.ndim = 2`` this is the usual diagonal, for
|
Chris@87
|
781 ``a.ndim > 2`` this is the set of indices to access ``a[i, i, ..., i]``
|
Chris@87
|
782 for ``i = [0..n-1]``.
|
Chris@87
|
783
|
Chris@87
|
784 Parameters
|
Chris@87
|
785 ----------
|
Chris@87
|
786 n : int
|
Chris@87
|
787 The size, along each dimension, of the arrays for which the returned
|
Chris@87
|
788 indices can be used.
|
Chris@87
|
789
|
Chris@87
|
790 ndim : int, optional
|
Chris@87
|
791 The number of dimensions.
|
Chris@87
|
792
|
Chris@87
|
793 See also
|
Chris@87
|
794 --------
|
Chris@87
|
795 diag_indices_from
|
Chris@87
|
796
|
Chris@87
|
797 Notes
|
Chris@87
|
798 -----
|
Chris@87
|
799 .. versionadded:: 1.4.0
|
Chris@87
|
800
|
Chris@87
|
801 Examples
|
Chris@87
|
802 --------
|
Chris@87
|
803 Create a set of indices to access the diagonal of a (4, 4) array:
|
Chris@87
|
804
|
Chris@87
|
805 >>> di = np.diag_indices(4)
|
Chris@87
|
806 >>> di
|
Chris@87
|
807 (array([0, 1, 2, 3]), array([0, 1, 2, 3]))
|
Chris@87
|
808 >>> a = np.arange(16).reshape(4, 4)
|
Chris@87
|
809 >>> a
|
Chris@87
|
810 array([[ 0, 1, 2, 3],
|
Chris@87
|
811 [ 4, 5, 6, 7],
|
Chris@87
|
812 [ 8, 9, 10, 11],
|
Chris@87
|
813 [12, 13, 14, 15]])
|
Chris@87
|
814 >>> a[di] = 100
|
Chris@87
|
815 >>> a
|
Chris@87
|
816 array([[100, 1, 2, 3],
|
Chris@87
|
817 [ 4, 100, 6, 7],
|
Chris@87
|
818 [ 8, 9, 100, 11],
|
Chris@87
|
819 [ 12, 13, 14, 100]])
|
Chris@87
|
820
|
Chris@87
|
821 Now, we create indices to manipulate a 3-D array:
|
Chris@87
|
822
|
Chris@87
|
823 >>> d3 = np.diag_indices(2, 3)
|
Chris@87
|
824 >>> d3
|
Chris@87
|
825 (array([0, 1]), array([0, 1]), array([0, 1]))
|
Chris@87
|
826
|
Chris@87
|
827 And use it to set the diagonal of an array of zeros to 1:
|
Chris@87
|
828
|
Chris@87
|
829 >>> a = np.zeros((2, 2, 2), dtype=np.int)
|
Chris@87
|
830 >>> a[d3] = 1
|
Chris@87
|
831 >>> a
|
Chris@87
|
832 array([[[1, 0],
|
Chris@87
|
833 [0, 0]],
|
Chris@87
|
834 [[0, 0],
|
Chris@87
|
835 [0, 1]]])
|
Chris@87
|
836
|
Chris@87
|
837 """
|
Chris@87
|
838 idx = arange(n)
|
Chris@87
|
839 return (idx,) * ndim
|
Chris@87
|
840
|
Chris@87
|
841
|
Chris@87
|
842 def diag_indices_from(arr):
|
Chris@87
|
843 """
|
Chris@87
|
844 Return the indices to access the main diagonal of an n-dimensional array.
|
Chris@87
|
845
|
Chris@87
|
846 See `diag_indices` for full details.
|
Chris@87
|
847
|
Chris@87
|
848 Parameters
|
Chris@87
|
849 ----------
|
Chris@87
|
850 arr : array, at least 2-D
|
Chris@87
|
851
|
Chris@87
|
852 See Also
|
Chris@87
|
853 --------
|
Chris@87
|
854 diag_indices
|
Chris@87
|
855
|
Chris@87
|
856 Notes
|
Chris@87
|
857 -----
|
Chris@87
|
858 .. versionadded:: 1.4.0
|
Chris@87
|
859
|
Chris@87
|
860 """
|
Chris@87
|
861
|
Chris@87
|
862 if not arr.ndim >= 2:
|
Chris@87
|
863 raise ValueError("input array must be at least 2-d")
|
Chris@87
|
864 # For more than d=2, the strided formula is only valid for arrays with
|
Chris@87
|
865 # all dimensions equal, so we check first.
|
Chris@87
|
866 if not alltrue(diff(arr.shape) == 0):
|
Chris@87
|
867 raise ValueError("All dimensions of input must be of equal length")
|
Chris@87
|
868
|
Chris@87
|
869 return diag_indices(arr.shape[0], arr.ndim)
|