Generate a n-dimensional array of coordinates in numpy

Question

Let's say I have a function f which can take coordinates as parameter and returns an integer (f(x) in this case). The coordinates can be multidimensional and are in the form of a list. My goal is to fill a numpy array with all values between two coordinates. I've tried to make a list of all possible indices and use it as input for the vectorized function.

Here is my code for 2 dimensional coordinates:

import itertools
import numpy


def index_array(lower_corner, upper_corner):
     x_range = range(lower_corner[0], upper_corner[0])
     y_range = range(lower_corner[1], upper_corner[1])
     return numpy.array(list(itertools.product(x_range, y_range)))


print(index_array([2, -2], [5, 3]))

This will return the index list like expected:

[[ 2 -2]
 [ 2 -1]
 [ 2  0]
 [ 2  1]
 [ 2  2]
 [ 3 -2]
 [ 3 -1]
 [ 3  0]
 [ 3  1]
 [ 3  2]
 [ 4 -2]
 [ 4 -1]
 [ 4  0]
 [ 4  1]
 [ 4  2]]

And here is my attempt for n dimensions:

import itertools
import numpy


def f(x):
    # dummy function
    return x + 5


def index_array(lower_corner, upper_corner):
    # returns all indices between two n-dimensional points
    range_list = []
    for n in range(len(lower_corner)):
        range_list.append(range(lower_corner[n], upper_corner[n]))
    return numpy.array(list(itertools.product(*range_list)))


lower_corner = numpy.array([2, -2])
upper_corner = numpy.array([5, 3])
indices = index_array(lower_corner, upper_corner)
vect_func = numpy.vectorize(f)
results = vect_func(indices)
print(results)

While this works it's quite slow and needs huge amounts of memory. Is it possible to write this in a more efficient way? I could think about using numpy.meshgrid but I don't know how I would use it.

Divakar · Accepted Answer · 2016-07-03 13:28:02Z

Indeed np.meshgrid would be one way to do it with some stacking, as shown below -

def ndim_grid(start,stop):
    # Set number of dimensions
    ndims = len(start)

    # List of ranges across all dimensions
    L = [np.arange(start[i],stop[i]) for i in range(ndims)]

    # Finally use meshgrid to form all combinations corresponding to all 
    # dimensions and stack them as M x ndims array
    return np.hstack((np.meshgrid(*L))).swapaxes(0,1).reshape(ndims,-1).T

Sample run

1) 2D Case :

In [97]: ndim_grid([2, -2],[5, 3])
Out[97]: 
array([[ 2, -2],
       [ 2, -1],
       [ 2,  0],
       [ 2,  1],
       [ 2,  2],
       [ 3, -2],
       [ 3, -1],
       [ 3,  0],
       [ 3,  1],
       [ 3,  2],
       [ 4, -2],
       [ 4, -1],
       [ 4,  0],
       [ 4,  1],
       [ 4,  2]])

2) 3D Case :

In [98]: ndim_grid([2, -2, 4],[5, 3, 6])
Out[98]: 
array([[ 2, -2,  4],
       [ 2, -2,  5],
       [ 2, -1,  4],
       [ 2, -1,  5],
       [ 2,  0,  4],
       [ 2,  0,  5],
       [ 2,  1,  4],
       [ 2,  1,  5],
       [ 2,  2,  4],
       [ 2,  2,  5],
       [ 3, -2,  4],
       [ 3, -2,  5],
       [ 3, -1,  4],
       [ 3, -1,  5],
       [ 3,  0,  4],
       [ 3,  0,  5],
       [ 3,  1,  4],
       [ 3,  1,  5],
       [ 3,  2,  4],
       [ 3,  2,  5],
       [ 4, -2,  4],
       [ 4, -2,  5],
       [ 4, -1,  4],
       [ 4, -1,  5],
       [ 4,  0,  4],
       [ 4,  0,  5],
       [ 4,  1,  4],
       [ 4,  1,  5],
       [ 4,  2,  4],
       [ 4,  2,  5]])

akuiper · Accepted Answer · 2016-07-03 13:39:22Z

3

Another option is to use the product from itertools, this also works if the corners are higher than 2D:

import itertools as it
lower_corner = [2, -2]
upper_corner = [5, 3]
[coord for coord in it.product(*[range(r[0], r[1]) for r in zip(lower_corner, upper_corner)])]

[(2, -2),
 (2, -1),
 (2, 0),
 (2, 1),
 (2, 2),
 (3, -2),
 (3, -1),
 (3, 0),
 (3, 1),
 (3, 2),
 (4, -2),
 (4, -1),
 (4, 0),
 (4, 1),
 (4, 2)]

answered Jul 3, 2016 at 13:39

akuiper

216k33 gold badges363 silver badges380 bronze badges

4 Comments

Gnarflord Over a year ago

Erm, I'm already using product from itertools . It's hidden in the return of the second function :D. But thanks for the more compact for loop!

Bakuriu Over a year ago

You could change range(r[0], r[1]) with range(*r). also you are using a list-comperhension that gets immediately unpacked. You could use a gen exp instead (i.e. *(range(*r) ... ) replacing [] with ()) and avoid creating that temporary list.

akuiper Over a year ago

@Bakuriu Works perfectly, thanks for pointing that out. A much more efficient and concise solution I believe.

aeolus Over a year ago

ty, this is such an elegant and easy solution for this problem. I'll add this: if one already has a list of vectors, just substitute that list into i.e. it.product(*my_vector_list) and ditch the even-spacing generator/list. This generalizes the solution or non-uniformly spaces vectors.

Collectives™ on Stack Overflow

Generate a n-dimensional array of coordinates in numpy

2 Answers 2

Comments

4 Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

2 Answers 2

Comments

4 Comments

Your Answer

Sign up or log in

Post as a guest

Related