[Numpy-discussion] numpy-izing a loop

[Robert Kern]

On Tue, Feb 10, 2009 at 10:19, Stéfan van der Walt <stefan at sun.ac.za> wrote:
> Hi Paul
>
> 2009/2/10 Paul Rudin <paul at rudin.co.uk>:
>>
>> I've just written this snippet of code:
>>
>> result = numpy.empty((dim, dim, dim), numpy.bool)
>> for x in xrange(dim):
>>    for y in xrange(dim):
>>        for z in xrange(dim):
>>            result[x, y, z] = ((xnear[y, z] < xfar[y, z]) and
>>                               (ynear[x, z] < yfar[x, z]) and
>>                               (znear[x, y] < zfar[x, y]))
>>
>>
>> Is there a way to get numpy to do the looping?
>
> You can try the following:
>
> x = np.arange(dim)
> y = np.arange(dim)
> z = np.arange(dim)
>
> result = (xnear[y,z] < xfar[y, z]) & (ynear[x, z] < yfar[x, z]) &
> (znear[x,y] < zfar[x,y])
>
> It broadcasts correctly, but you'll have to verify the results to be sure.

It doesn't broadcast correctly. You will get result.shape == (dim,)
but he wants (dim,dim,dim). You can get this with the following:

  x = np.arange(dim)[:,np.newaxis,np.newaxis]
  y = np.arange(dim)[np.newaxis,:,np.newaxis]
  z = np.arange(dim)[np.newaxis,np.newaxis,:]

-- 
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless
enigma that is made terrible by our own mad attempt to interpret it as
though it had an underlying truth."
  -- Umberto Eco