Skip to content

Instantly share code, notes, and snippets.

@zarch

zarch/resample_array.py

Created March 26, 2014 13:18
Show Gist options
  • Save zarch/9782938 to your computer and use it in GitHub Desktop.
Save zarch/9782938 to your computer and use it in GitHub Desktop.
Download ZIP
resample a numpy array
Raw
resample_array.py
# -*- coding: utf-8 -*-
"""
Created on Tue Mar 25 22:55:20 2014
@author: pietro
"""
import itertools
import numpy as np
def cartesian(arrays, out=None):
"""
Generate a cartesian product of input arrays.
Parameters
----------
arrays : list of array-like
1-D arrays to form the cartesian product of.
out : ndarray
Array to place the cartesian product in.
Returns
-------
out : ndarray
2-D array of shape (M, len(arrays)) containing cartesian products
formed of input arrays.
Examples
--------
::
>>> cartesian(([1, 2, 3], [4, 5], [6, 7]))
array([[1, 4, 6],
[1, 4, 7],
[1, 5, 6],
[1, 5, 7],
[2, 4, 6],
[2, 4, 7],
[2, 5, 6],
[2, 5, 7],
[3, 4, 6],
[3, 4, 7],
[3, 5, 6],
[3, 5, 7]])
Reference
---------
http://stackoverflow.com/a/1235363
"""
arrays = [np.asarray(x) for x in arrays]
dtype = arrays[0].dtype
n = np.prod([x.size for x in arrays])
if out is None:
out = np.zeros([n, len(arrays)], dtype=dtype)
m = n / arrays[0].size
out[:, 0] = np.repeat(arrays[0], m)
if arrays[1:]:
cartesian(arrays[1:], out=out[0:m, 1:])
for j in range(1, arrays[0].size):
out[j*m:(j+1)*m, 1:] = out[0:m, 1:]
return out
def split_in_chunk(iterable, lenght=10):
"""Split a iterable object in chunks.
Parameters
----------
iterable : iterable
A iterable object
lenght : integer
Max lenght of each chunk.
Returns
-------
out : iterator
Iterator of chunks.
Examples
--------
::
>>> for chunk in split_in_chunk(range(25)):
... print(chunk)
(0, 1, 2, 3, 4, 5, 6, 7, 8, 9)
(10, 11, 12, 13, 14, 15, 16, 17, 18, 19)
(20, 21, 22, 23, 24)
>>> for chunk in split_in_chunk(range(25), 3):
... print(chunk)
(0, 1, 2)
(3, 4, 5)
(6, 7, 8)
(9, 10, 11)
(12, 13, 14)
(15, 16, 17)
(18, 19, 20)
(21, 22, 23)
(24,)
"""
it = iter(iterable)
while True:
chunk = tuple(itertools.islice(it, lenght))
if not chunk:
return
yield chunk
def get_chunks(ashape, tshape):
"""Return a tuple with a list of arrays with chunks of indexes.
Parameters
----------
ashape : shape of the array
A tuple with the number of rows and the number of columns
tshape : shape of the tile
A tuple with the number of rows and the number of columns
Returns
-------
out : tuple
A tuple wiht two list of arrays with the indexes
Examples
--------
::
>>> ix, iy = get_chunks((5, 5), (2, 2))
>>> ix
[array([0, 1]), array([2, 3]), array([4])]
>>> iy
[array([0, 1]), array([2, 3]), array([4])]
"""
nrows, ncols = ashape
trows, tcols = tshape
return ([np.array(chnk) for chnk in split_in_chunk(range(nrows), trows)],
[np.array(chnk) for chnk in split_in_chunk(range(ncols), tcols)])
def resample_array(array, tile_shape, default=0):
"""Split a 2D array into tiles, returning a 3D array.
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| 1, 1 | 2, 2 | 3 |
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| 1, 1 | 2, 2 | 3 |
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| | | |
+------+------+------+
Parameters
----------
array : ndarray
2-D array
tile_shape : a tuple
A tuple with the number of rows and columns of the tile.
Returns
-------
out : ndarray
3D ndarray
Examples
--------
::
>>> array = np.array([[1, 1, 2, 2, 3, ], ] * 5)
>>> array
array([[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3]])
>>> resample_array(array, (2, 2))
array([[[1, 1, 1, 1],
[2, 2, 2, 2],
[3, 3, 0, 0]],
<BLANKLINE>
[[1, 1, 1, 1],
[2, 2, 2, 2],
[3, 3, 0, 0]],
<BLANKLINE>
[[1, 1, 0, 0],
[2, 2, 0, 0],
[3, 0, 0, 0]]])
>>> np.sum(resample_array(array, (2, 2)), axis=2)
array([[4, 8, 6],
[4, 8, 6],
[2, 4, 3]])
>>> np.mean(resample_array(array, (2, 2)), axis=2)
array([[ 1. , 2. , 1.5 ],
[ 1. , 2. , 1.5 ],
[ 0.5 , 1. , 0.75]])
"""
rchks, cchks = get_chunks(array.shape, tile_shape)
rows = len(rchks)
cols = len(cchks)
result = np.zeros((rows, cols, tile_shape[0] * tile_shape[1]),
dtype=array.dtype)
if default != 0:
result[:] = default
for i in range(rows):
for j in range(cols):
ix, iy = cartesian([rchks[i], cchks[j]]).T
result[i, j, :len(ix)] = array[ix, iy]
return result
def resample(array, tile_shape, function, default=0):
"""Split a 2D array into tiles, returning a 3D array.
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| 1, 1 | 2, 2 | 3 |
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| 1, 1 | 2, 2 | 3 |
+------+------+------+
| 1, 1 | 2, 2 | 3 |
| | | |
+------+------+------+
Parameters
----------
array : ndarray
2-D array
tile_shape : a tuple
A tuple with the number of rows and columns of the tile.
Returns
-------
out : ndarray
3D ndarray
Examples
--------
::
>>> array = np.array([[1, 1, 2, 2, 3, ], ] * 5)
>>> array
array([[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3],
[1, 1, 2, 2, 3]])
>>> resample(array, (2, 2), np.sum)
array([[4, 8, 6],
[4, 8, 6],
[2, 4, 3]])
>>> resample(array, (2, 2), np.mean)
array([[1, 2, 3],
[1, 2, 3],
[1, 2, 3]])
"""
rchks, cchks = get_chunks(array.shape, tile_shape)
rows = len(rchks)
cols = len(cchks)
result = np.zeros((rows, cols), dtype=array.dtype)
if default != 0:
result[:] = default
for i in range(rows):
for j in range(cols):
ix, iy = cartesian([rchks[i], cchks[j]]).T
result[i, j] = function(array[ix, iy])
return result
if __name__ == "__main__":
import doctest
doctest.testmod()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment