[Numpy-discussion] Updating SciPy from 0.9.0 to 0.10.1 triggers undesired behavior in NumPy 1.6.2

[Wolfgang Draxinger]

Hi,

first my apologies for crossposting this to two maillists, but as both
projects are affected I think this is in order.

Like the subject says I encountered some undesired behavior in the
interaction of SciPy with NumPy.

Using the "old" SciPy version 0.9.0 everything works fine and smooth.
But upgrading to SciPy-0.10.1 triggers some ValueError in
numpy.histogram2d of NumPy-1.6.2 when executing one of my programs.

I developed a numerical evaluation system for our detectors here. One of
the key operations is determining the distribution of some 2-dimensional
variable space based on the values found in the image delivered by the
detector, where each pixel has associated values for the target
variables. This goes something like the following

   ABdist, binsA, binsB = numpy.histogram2d(
       B_yz.ravel(), 
       A_yz.ravel(), 
       [B_bins, A_bins],
       weights=image.ravel() )

The bins parameter can be either [int, int] or [array, array], that
makes no difference in the outcome.

The mappings A_yz and B_yz are created using
scipy.interpolate.griddata. We have a list of pairs of pairs which are
determined by measurement. Basically in the calibration step we vary
variables A,B and store at which Y,Z we get the corresponding signal.
So essentially this is a (A,B) -> (Y,Z) mapping. In the region of
interest is has a bijective subset that's also smooth. However the
original data also contains areas where (Y,Z) has no corresponding
(A,B) or where multiple (A,B) map to the same (Y,Z); like said, those
lie outside the RoI. For our measurements we need to reverse this
process, i.e. we want to do (Y,Z) -> (A,B).

So I use griddata to evaluate a discrete reversal for this
mapping, of the same dimensions that the to be evaluated image has:

   gry, grz = numpy.mgrid[self.y_lim[0]:self.y_lim[1]:self.y_res*1j,
                          self.z_lim[0]:self.z_lim[1]:self.z_res*1j]
   # for whatever reason I have to do the following
   # assigning to evalgrid directly breaks the program.
   evalgrid = (gry, grz)

   points = (Y.ravel(), Z.ravel())

   def gd(a):
       return scinp.griddata(
           points, 
           a.ravel(), 
           evalgrid, 
           method='cubic' )

   A_yz = gd(A)
   B_yz = gd(B)

where A,B,Y,Z have the same dimensions and are the ordered lists/arrays
of the scalar values of the two sets mapped between. As you can see,
this approach does also involve the elements of the sets, which are not
mapped bijectively. As lying outside the convex boundary or not being
properly interpolatable they should receive the fill value.

As long as I stay with SciPy-0.9.0 everything works fine. However after
upgrading to SciPy-0.10.1 the histogram2d step fails with a ValueError.
The version of NumPy is 1.6.2 for both cases.

/usr/lib64/python2.7/site-packages/numpy/ma/core.py:772: RuntimeWarning: invalid value encountered in absolute
  return umath.absolute(a) * self.tolerance >= umath.absolute(b)
Traceback (most recent call last):
  File "./ephi.py", line 71, in <module>
    ABdist, binsA, binsB = numpy.histogram2d(B_yz.ravel(), A_yz.ravel(), [B_bins, A_bins], weights=image.ravel())
  File "/usr/lib64/python2.7/site-packages/numpy/lib/twodim_base.py", line 615, in histogram2d
    hist, edges = histogramdd([x,y], bins, range, normed, weights)
  File "/usr/lib64/python2.7/site-packages/numpy/lib/function_base.py", line 357, in histogramdd
    decimal = int(-log10(mindiff)) + 6
ValueError: cannot convert float NaN to integer

Any ideas on this?


Regards,

Wolfgang
-- 
Fakultät für Physik, LMU München
Beschleunigerlabor/MLL
Am Coulombwall 6
85748 Garching
Deutschland

Tel: +49-89-289-14286
Fax: +49-89-289-14280