floating point bug

[arushanova]

@jwaterfield I've checked the shift with numbers, that caused an error previously, now its ok!

I've checked shifting the spectrum completely outside the boundaries, i.e. there were no events left at all. The method plot_root.plot_projection just returned an empty canvas, while the method plot.plot_projection was returning this error message:

fig_1 = plot.plot_projection(very_new, "energy_mc", fig_num=1)

File "/.data/snoplus/evelina/virtual_echidna/echidna/echidna/output/plot.py", line 75, in plot_projection axis.hist(x, bins, weights=data, histtype="step") File "/.data/snoplus/evelina/virtual_echidna/lib/python2.7/site-packages/matplotlib/axes/_axes.py", line 5874, in hist ymin = max(ymin*0.9, minimum) UnboundLocalError: local variable 'ymin' referenced before assignment

I believe we want to prevent shiftind the spectrum outside the boundaries and quit the job with sufficient comment.

[ashleyrback]

I think this just needs to have something in the _shift and _shift_by_bin methods that raises an error if par._low + self._shift >= par._high or par._high + self._shift <= par._low.

[arushanova]

@jwaterfield I have similar error message if I use scale. So here is what I do: scale_class = scale.Scale() scale_class.set_scale_factor(0.9) scaled_spectrum = scale_class.scale(smeared_spectrum, "energy_mc") print scaled_spectrum._data plot.plot_projection(scaled_spectrum, "energy_mc")

For the data print out I have a list of nans: ...nan nan nan nan nan nan nan nan]

And here is the error message: File "/.data/snoplus/evelina/virtual_echidna/echidna/echidna/output/plot.py", line 75, in plot_projection axis.hist(x, bins, weights=data, histtype="step") File "/.data/snoplus/evelina/virtual_echidna/lib/python2.7/site-packages/matplotlib/axes/_axes.py", line 5874, in hist ymin = max(ymin*0.9, minimum) UnboundLocalError: local variable 'ymin' referenced before assignment

Since the error messages are similar, could you also update the scaling in this pull request?

[ashleyrback]

@jwaterfield, how are things going with this PR? Do you want me to make the changes (suggested above) to the Shift class, which you can them merge in. That way you only have to concentrate on the Scale class.

[arushanova]

@jwaterfield @ashleyrback

I've tested scale and shift. So,

1. Shift

When the spectrum is not "touching zero bin" - i works fine. The problem appears, when the spectrum starts at zero. In this case shift_by_bin works, but not a normal fift. I have 500 bins and 500 is max value, I was shifting in 0.1. The error message that I get:

Traceback (most recent call last): File "to_test_bugs.py", line 28, in plot.plot_projection(shifted_spectrum, "energy_mc") File "/.data/snoplus/evelina/virtual_echidna/echidna/echidna/output/plot.py", line 75, in plot_projection axis.hist(x, bins, weights=data, histtype="step") File "/.data/snoplus/evelina/virtual_echidna/lib/python2.7/site-packages/matplotlib/axes/_axes.py", line 5874, in hist ymin = max(ymin*0.9, minimum) UnboundLocalError: local variable 'ymin' referenced before assignment

The data has a list full of nans: nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan]

1. Scale scale_class = scale.Scale() scale_class.set_scale_factor( value) If I scale spectrum to a value greater than 1, than no matter where the spectrum is - it works. But if a value is less than 1, here is an error message:

Traceback (most recent call last): File "to_test_bugs.py", line 19, in spectrum = scale_class.scale(spectrum, "energy_mc") File "/.data/snoplus/evelina/virtual_echidna/echidna/echidna/core/scale.py", line 66, in scale old_bin1 = par.get_bin(x/sf) File "/.data/snoplus/evelina/virtual_echidna/echidna/echidna/core/spectra.py", line 150, in get_bin % (x, self._high)) ValueError: 10.0055555556 is above parameter upper bound 10.0

[jwaterfield]

@arushanova Would you be able to email this hdf5 and I'll try to fix it

[arushanova]

@jwaterfield I've just sent both files with I use for testing. I also use gaus, because it's easy to create it at any point along the x axis to test the "border" effect.

[arushanova]

@jwaterfield @ashleyrback @EdLeming

About the scaling. I ran debugger. So basically when you're scaling down and do x/sf - you're going outside the boundaries.

/.data/snoplus/evelina/virtual_echidna/echidna/echidna/core/scale.py(68)scale() -> print x, sf (Pdb) n -9.995 0.9 /.data/snoplus/evelina/virtual_echidna/echidna/echidna/core/scale.py(69)scale() -> old_bin1 = par.get_bin(x/sf) (Pdb) n ValueError: ValueErr... -10.0',) /.data/snoplus/evelina/virtual_echidna/echidna/echidna/core/scale.py(69)scale() -> old_bin1 = par.get_bin(x/sf) (Pdb) n --Return--

Also, I've changed the unit test, originaly you had scale =1.1, I changed it to 0.1 and it failed.

[arushanova]

@jwaterfield @EdLeming

Ashley and I are investigating this bug

[ashleyrback]

@jwaterfield, the latest version works fine for sf ~1, but @arushanova still sees problems for a scale factor of 0.01, for example. Whilst I agree, that in most cases we would not use a sf of 0.01, it has highlighted some parts of the code which concern me, in terms of the accuracy of scale convolutions that are in the ranges we will want to use.

The main problem happens in this line, because the scaled_spec.sum() == 0. Then when you set scaled_spec._num_decays equal to this value, subsequent scaling results in nan for all values - as you end up with a zero-division.

We traced this problem back to this line, where old_data1 evaluates as [ 0. ], as does old_data2 below, for the first bin. The slices it is taking are [49:50] and [50:51]. Evelina can perhaps paste the full _data array below, so you can see a bit better. But the real problem is that the for subsequent bins, the slices will move to the right, but all the events around bins 0-10 never get included. Meaning that at the end, the scaled spectrum contains no events.

As I siad, you don't see this problem for scaling factors ~1, since the slice never strays too far from the current bin.

I'll highlight, in the code some other concerns this has highlighted.

[arushanova]

This is an array Ashley was talking about: (Pdb) print spectrum._data [ 0. 9. 25. 62. 267. 906. 420. 398. 476. 571. 542. 603.

1. 656. 625. 591. 512. 429. 401. 320. 230. 176. 111. 71.
   657. 42. 17. 10. 4. 2. 0. 0. 0. 0. 0. 0.
   658. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0.
   659. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
   660. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
   661. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
   662. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
   663. 0. 0. 0.]

[ashleyrback]

@arushanova can you test the updates that James pushed here, with your script that was breaking before, to see if this fixes it please?