How to convert the output of meshgrid to the corresponding array of points?

Question

I want to create a list of points that would correspond to a grid. So if I want to create a grid of the region from (0, 0) to (1, 1), it would contain the points (0, 0), (0, 1), (1, 0) and (1, 0).

I know that that this can be done with the following code:

g = np.meshgrid([0,1],[0,1]) np.append(g[0].reshape(-1,1),g[1].reshape(-1,1),axis=1) 

Yielding the result:

array([[0, 0],        [1, 0],        [0, 1],        [1, 1]]) 

My question is twofold:

1. Is there a better way of doing this?
2. Is there a way of generalizing this to higher dimensions?

Answer 1

I just noticed that the documentation in numpy provides an even faster way to do this:

X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:100j] positions = np.vstack([X.ravel(), Y.ravel()]) 

This can easily be generalized to more dimensions using the linked meshgrid2 function and mapping 'ravel' to the resulting grid.

g = meshgrid2(x, y, z) positions = np.vstack(map(np.ravel, g)) 

The result is about 35 times faster than the zip method for a 3D array with 1000 ticks on each axis.

Source: http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.gaussian_kde.html#scipy.stats.gaussian_kde

To compare the two methods consider the following sections of code:

Create the proverbial tick marks that will help to create the grid.

In [23]: import numpy as np  In [34]: from numpy import asarray  In [35]: x = np.random.rand(100,1)  In [36]: y = np.random.rand(100,1)  In [37]: z = np.random.rand(100,1) 

Define the function that mgilson linked to for the meshgrid:

In [38]: def meshgrid2(*arrs):    ....:     arrs = tuple(reversed(arrs))    ....:     lens = map(len, arrs)    ....:     dim = len(arrs)    ....:     sz = 1    ....:     for s in lens:    ....:        sz *= s    ....:     ans = []    ....:     for i, arr in enumerate(arrs):    ....:         slc = [1]*dim    ....:         slc[i] = lens[i]    ....:         arr2 = asarray(arr).reshape(slc)    ....:         for j, sz in enumerate(lens):    ....:             if j != i:    ....:                 arr2 = arr2.repeat(sz, axis=j)    ....:         ans.append(arr2)    ....:     return tuple(ans) 

Create the grid and time the two functions.

In [39]: g = meshgrid2(x, y, z)  In [40]: %timeit pos = np.vstack(map(np.ravel, g)).T 100 loops, best of 3: 7.26 ms per loop  In [41]: %timeit zip(*(x.flat for x in g)) 1 loops, best of 3: 264 ms per loop