[Numpy-discussion] Is this a bug in numpy.ma.reduce?
David Goldsmith
d.l.goldsmith at gmail.com
Mon Mar 8 14:25:03 EST 2010
On Mon, Mar 8, 2010 at 10:17 AM, David Goldsmith <d.l.goldsmith at gmail.com>wrote:
> On Mon, Mar 8, 2010 at 6:52 AM, Bruce Southey <bsouthey at gmail.com> wrote:
>
>> On 03/08/2010 01:30 AM, David Goldsmith wrote:
>>
>> On Sun, Mar 7, 2010 at 4:41 AM, Friedrich Romstedt <
>> friedrichromstedt at gmail.com> wrote:
>>
>>> I would like to stress the fact that imo this is maybe not ticket and not
>>> a bug.
>>>
>>> The issue arises when calling a.max() or similar of empty arrays a, i.e.,
>>> with:
>>>
>>> >>> 0 in a.shape
>>> True
>>>
>>> Opposed to the .prod() of an empty array, such a .max() or .min()
>>> cannot be defined, because the set is empty. So it's fully correct to
>>> let such calls fail. Just the failure is a bit deep in numpy, and
>>> only the traceback gives some hint what went wrong.
>>>
>>> I posted something similar also on the matplotlib-users list, sorry
>>> for cross-posting thus.
>>>
>>
>> Any suggestions, then, how to go about figuring out what's happening in my
>> code that's causing this "feature" to manifest itself?
>>
>> DG
>>
>> Perhaps providing the code with specific versions of Python, numpy etc.
>> would help.
>>
>> I would guess that aquarius_test.py has not correctly setup the necessary
>> inputs (or has invalid inputs) required by matplotlib (which I have no
>> knowledge about). Really you have to find if the _A in cmp.py used by
>> 'self.norm.autoscale_None(self._A)' is valid. You may be missing a valid
>> initialization step because the TypeError exception in autoscale_None ('You
>> must first set_array for mappable') implies something need to be done first.
>>
>>
>> Bruce
>>
>
> Python 2.5.4, Numpy 1.4.0, Matplotlib 0.99.0, Windows 32bit Vista Home
> Premium SP2
>
> # Code copyright 2010 by David Goldsmith
> # Comments and unnecessaries edited for brevity
> import numpy as N
> import matplotlib as MPL
> from matplotlib import pylab
> from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
> from matplotlib.figure import Figure
> import matplotlib.cm as cm
>
> J = complex(0,1); tol = 1e-6; maxiter = 20;
> roots = (-2 + J/2, -1 + J, -0.5 + J/2, 0.5 + J, 1 + J/2, 2 + J, 2.5 +
> J/2,
> -2 - J, -1 - J/2, -0.5 - J, 0.5 - J/2, 1 - J, 2 - J/2, 2.5 -
> J)
>
> def ffp(z):
> w, wp = (0J, 0J)
> for root in roots:
> z0 = z - root
> w += N.sin(1/z0)
> wp -= N.cos(1/z0)/(z0*z0)
> return (w, wp)
>
> def iter(z):
> w, wp = ffp(z)
> return z - w/wp
>
> def find_root(z0):#, k, j):
> count = 0
> z1 = iter(z0)
> if N.isnan(z1):
> return N.complex64(N.inf)
> while (N.abs(z1 - z0) > tol) and \
> (count < maxiter):
> count += 1
> z0 = z1
> z1 = iter(z0)
> if N.abs(z1 - z0) > tol:
> result = 0
> else:
> result = z1
> return N.complex64(result)
>
> w, h, DPI = (3.2, 2.0, 100)
> fig = Figure(figsize=(w, h),
> dpi=DPI,
> frameon=False)
> ax = fig.add_subplot(1,1,1)
> canvas = FigureCanvas(fig)
> nx, xmin, xmax = (int(w*DPI), -0.5, 0.5)
> ny, ymin, ymax = (int(h*DPI), 0.6, 1.2)
>
> X, xincr = N.linspace(xmin,xmax,nx,retstep=True)
> Y, yincr = N.linspace(ymin,ymax,ny,retstep=True)
> W = N.zeros((ny,nx), dtype=N.complex64)
>
> for j in N.arange(nx):
> if not (j%100): # report progress
> print j
> for k in N.arange(ny):
> x, y = (X[j], Y[k])
> z0 = x + J*y
> W[k,j] = find_root(z0)#,k,j)
>
> print N.argwhere(N.logical_not(N.isfinite(W.real)))
> print N.argwhere(N.logical_not(N.isfinite(W.imag)))
> W = W.T
> argW = N.angle(W)
> print N.argwhere(N.logical_not(N.isfinite(argW)))
> cms = ("Blues",)# "Blues_r", "cool", "cool_r",
>
> def all_ticks_off(ax):
> ax.xaxis.set_major_locator(pylab.NullLocator())
> ax.yaxis.set_major_locator(pylab.NullLocator())
>
> for cmap_name in cms:
> all_ticks_off(ax)
> ax.hold(True)
> for i in range(4):
> for j in range(4):
> part2plot = argW[j*ny/4:(j+1)*ny/4, i*nx/4:(i+1)*nx/4]
> if N.any(N.logical_not(N.isfinite(part2plot))):
> print i, j,
> print N.argwhere(N.logical_not(N.isfinite(part2plot)))
> extent = (i*nx/4, (i+1)*nx/4, (j+1)*ny/4, j*ny/4)
>
> ax.imshow(part2plot, cmap_name, extent = extent)
> ax.set_xlim(0, nx)
> ax.set_ylim(0, ny)
> canvas.print_figure('../../Data-Figures/Zodiac/Aquarius/'+ cmap_name +
> 'Aquarius_test.png', dpi=DPI)
> # End Aquarius_test.png
>
> DG
>
Oh, and here's "fresh" output (i.e., I just reran it to confirm that I'm
still having the problem).
0
100
200
300
[[133 319]]
[]
[]
Traceback (most recent call last):
File
"C:\Users\Fermat\Documents\Fractals\Python\Source\Zodiac\aquarius_test.py",
line 108, in <module>
ax.imshow(part2plot, cmap_name, extent = extent)
File "C:\Python254\lib\site-packages\matplotlib\axes.py", line 6261, in
imshow
im.autoscale_None()
File "C:\Python254\lib\site-packages\matplotlib\cm.py", line 236, in
autoscale_None
self.norm.autoscale_None(self._A)
File "C:\Python254\lib\site-packages\matplotlib\colors.py", line 792, in
autoscale_None
if self.vmin is None: self.vmin = ma.minimum(A)
File "C:\Python254\Lib\site-packages\numpy\ma\core.py", line 5555, in
__call__
return self.reduce(a)
File "C:\Python254\Lib\site-packages\numpy\ma\core.py", line 5570, in reduce
t = self.ufunc.reduce(target, **kargs)
ValueError: zero-size array to ufunc.reduce without identity
Script terminated.
Notes:
0) Despite
N.argwhere(N.logical_not(N.isfinite(W.real)))
being non-empty,
N.argwhere(N.logical_not(N.isfinite(N.angle(W.T))))
apparently is empty, so I suppose whatever is behind this may also be the
source of my problem, but:
1) Nothing is ever printed as a result of this attempt to catch a non-finite
problem:
if N.any(N.logical_not(N.isfinite(part2plot))):
print i, j,
print N.argwhere(N.logical_not(N.isfinite(part2plot)))
which occurs for all part2plot's before they're passed to imshow, so I
conclude that the problem is something other than non-finites?
DG
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20100308/473500be/attachment.html>
More information about the NumPy-Discussion
mailing list