Split a 3D numpy array into 3D blocks
Question: I would like to split a 3D numpy array into 3D blocks in a 'pythonic' way. I
am working with image sequences that are somewhat large arrays
(1000X1200X1600), so I need to split them into pieces to do my processing.

I have written functions to do this, but I am wondering if there is a native
numpy way to accomplish this - numpy.split does not seem to do what I want for
3D arrays (but perhaps I don't understand its functionality)

To be clear: the code below accomplishes my task, but I am seeking a faster
way to do it.

    
    
    def make_blocks(x,t):
    #x should be a yXmXn matrix, and t should even divides m,n
    #returns a list of 3D blocks of size yXtXt 
        down =  range(0,x.shape[1],t)
        across = range(0,x.shape[2],t)
        reshaped = []
        for d in down:
            for a in across:
                reshaped.append(x[:,d:d+t,a:a+t])
        return reshaped
    
    def unmake_blocks(x,d,m,n):
    #this takes a list of matrix blocks of size dXd that is m*n/d^2 long 
    #returns a 2D array of size mXn
        rows = []
        for i in range(0,int(m/d)):
            rows.append(np.hstack(x[i*int(n/d):(i+1)*int(n/d)]))
        return np.vstack(rows)
    


Answer: Here are vectorized versions of those loopy implementations using a
combination of permuting dims with
[`np.transpose`](http://docs.scipy.org/doc/numpy/reference/generated/numpy.transpose.html)
and
[`reshaping`](http://docs.scipy.org/doc/numpy/reference/generated/numpy.reshape.html)
-

    
    
    def make_blocks_vectorized(x,d):
        p,m,n = x.shape
        return x.reshape(-1,m//d,d,n//d,d).transpose(1,3,0,2,4).reshape(-1,p,d,d)
    
    def unmake_blocks_vectorized(x,d,m,n):    
        return np.concatenate(x).reshape(m//d,n//d,d,d).transpose(0,2,1,3).reshape(m,n)
    

* * *

Sample run for `make_blocks` -

    
    
    In [120]: x = np.random.randint(0,9,(2,4,4))
    
    In [121]: make_blocks(x,2)
    Out[121]: 
    [array([[[4, 7],
             [8, 3]],
    
            [[0, 5],
             [3, 2]]]), array([[[5, 7],
             [4, 0]],
    
            [[7, 3],
             [5, 7]]]), ... and so on.
    
    In [122]: make_blocks_vectorized(x,2)
    Out[122]: 
    array([[[[4, 7],
             [8, 3]],
    
            [[0, 5],
             [3, 2]]],
    
    
           [[[5, 7],
             [4, 0]],
    
            [[7, 3],
             [5, 7]]],  ... and so on.
    

Sample run for `unmake_blocks` -

    
    
    In [135]: A = [np.random.randint(0,9,(3,3)) for i in range(6)]
    
    In [136]: d = 3
    
    In [137]: m,n = 6,9
    
    In [138]: unmake_blocks(A,d,m,n)
    Out[138]: 
    array([[6, 6, 7, 8, 6, 4, 5, 4, 8],
           [8, 8, 3, 2, 7, 6, 8, 5, 1],
           [5, 2, 2, 7, 1, 2, 3, 1, 5],
           [6, 7, 8, 2, 2, 1, 6, 8, 4],
           [8, 3, 0, 4, 4, 8, 8, 6, 3],
           [5, 5, 4, 8, 5, 2, 2, 2, 3]])
    
    In [139]: unmake_blocks_vectorized(A,d,m,n)
    Out[139]: 
    array([[6, 6, 7, 8, 6, 4, 5, 4, 8],
           [8, 8, 3, 2, 7, 6, 8, 5, 1],
           [5, 2, 2, 7, 1, 2, 3, 1, 5],
           [6, 7, 8, 2, 2, 1, 6, 8, 4],
           [8, 3, 0, 4, 4, 8, 8, 6, 3],
           [5, 5, 4, 8, 5, 2, 2, 2, 3]])
    

* * *

Alternative to `make_blocks` with [`view_as_blocks`](http://scikit-
image.org/docs/dev/api/skimage.util.html#skimage.util.view_as_blocks) -

    
    
    from skimage.util.shape import view_as_blocks
    
    def make_blocks_vectorized_v2(x,d):
        return view_as_blocks(x,(x.shape[0],d,d))
    

* * *

**Runtime test**

1) `make_blocks` with original and `view_as_blocks` based approaches -

    
    
    In [213]: x = np.random.randint(0,9,(100,160,120)) # scaled down by 10
    
    In [214]: %timeit make_blocks(x,10)
    1000 loops, best of 3: 198 µs per loop
    
    In [215]: %timeit view_as_blocks(x,(x.shape[0],10,10))
    10000 loops, best of 3: 85.4 µs per loop
    

2) `unmake_blocks` with original and `transpose+reshape` based approaches -

    
    
    In [237]: A = [np.random.randint(0,9,(10,10)) for i in range(600)]
    
    In [238]: d = 10
    
    In [239]: m,n = 10*20,10*30
    
    In [240]: %timeit unmake_blocks(A,d,m,n)
    100 loops, best of 3: 2.03 ms per loop
    
    In [241]: %timeit unmake_blocks_vectorized(A,d,m,n)
    1000 loops, best of 3: 511 µs per loop
    

