[Rivet] Two-dimensional histograms

[Andy Buckley]

Hi James,

Weird! The condition that produces that error is the following:

if (b_acc.effNumEntries() > b_tot.effNumEntries() || b_acc.sumW() >
b_tot.sumW()) // -> throw

The effNumEntries is sum(W)**2 / sum(W**2) so that if the weights are
uniform (but not necessarily unity) it is equal to the number of fills.
If the weight distribution isn't uniform, it will differ and is more
robust than the raw number of entries. Note that the condition above
doesn't throw in the case of 100% efficiency, i.e. exact same fill stats
in accepted & total histos.

I don't see an error in the code, and we do handle the case (which you
have in your last few examples) where the total histo is unfilled --
although I would like to do something nicer than always set the result
to zero in that undefined case. Any ideas what could be going on? Are
your weights uniform?

Thanks,
Andy


On 23/09/14 20:16, James Robinson wrote:
> Dear Andy,
> 
> I'm getting an error from Rivet complaining about an invalid efficiency
> calculation:
> 
> Rivet_i              INFO Rivet_i finalizing
> Rivet.Analysis.Handler: INFO  Finalising analyses
> Rivet_i             FATAL  Standard std::exception is caught
> Rivet_i             ERROR Attempt to calculate an efficiency when the
> numerator is not a subset of the denominator
> 
> However, if I get the code to print out the individual histograms,
> bin-by-bin, I can't see any bin where the numerator is larger than the
> denominator:
> 
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  There are 263 / 393 entri
> ​es in the histograms
> ​​
> ​ ​
> 
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 1.1726e+17, total:
> 1.23204e+17
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 1.01373e+17,
> total: 1.23149e+17
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 6.85239e+16,
> total: 1.2764e+17
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 8.20248e+16,
> total: 1.4139e+17
> ​
> 
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 1.65158e+16,
> total: 4.2562e+16
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 7.17599e+15,
> total: 1.69049e+16
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 7.2618e+14, total:
> 2.16647e+16
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 9.2832e+15, total:
> 2.75894e+16
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 2.69685e+14,
> total: 3.10438e+14
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 9.55894e+14,
> total: 9.55894e+14
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 0, total: 1.66588e+15
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 0, total: 0
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 0, total: 0
> Rivet.Analysis.ATLAS_STDM_2012_17: INFO  dy accepted: 0, total: 0
> 
> Any idea what I'm doing wrong?
> 
> Thanks,
> James
> 
> On 17 September 2014 11:13, Andy Buckley <andy.buckley at cern.ch
> <mailto:andy.buckley at cern.ch>> wrote:
> 
>     On 17/09/14 10:38, James Robinson wrote:
>     >
>     >
>     > On 17 September 2014 11:26, Andy Buckley <andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>
>     > <mailto:andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>>> wrote:
>     >
>     >     On 16/09/14 22:53, James Robinson wrote:
>     >     > Dear Andy and David,
>     >     >
>     >     > On 16 September 2014 20:58, Andy Buckley <andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>
>     <mailto:andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>>
>     >     > <mailto:andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>
>     <mailto:andy.buckley at cern.ch <mailto:andy.buckley at cern.ch>>>> wrote:
>     >     >
>     >     >     In fact, there is even an efficiency(pass, tot) function that takes 1D
>     >     >     histos and returns a Scatter2D with appropriate binomial statistics
>     >     >     treatment!
>     >     >
>     >     >     Andy
>     >     >
>     >     >
>     >     > ​I think this is exactly what I want! The binomial error treatment is
>     >     > important for this analysis as the efficiencies tend towards 1.0 in some
>     >     > parts of phase space.​ What is the appropriate class called?
>     >
>     >     It's a function rather than a bound method (so call it as above):
>     >
>     >     https://yoda.hepforge.org/trac/browser/include/YODA/Histo1D.h#L426
>     >
>     >     We don't yet have the equivalent for 2D histos, but it should be trivial
>     >     to put together from the 1D version. I'll add that to my TODO list, but
>     >     maybe someone else fancies doing it? ;-)
>     >
>     >
>     > ​Ah OK. That makes sense. Presumably I still need to book it somehow to
>     > get it written out? Is there another method to register existing YODA
>     > objects? Or can the bookScatter2D() function accept an existing
>     > Scatter2D as it's argument (instead of a bin descriptor)?
> 
>     Good points -- I should provide helper methods for doing this in Rivet,
>     with automatic registration of the output Scatter, like we have for
>     YODA's divide(h_a,h_b) -> Rivet's divide(h_a,h_b,s_out). I'll add those
>     to the 2.2.0 release candidate.
> 
>     For now I think you can use addAnalysisObject(aoptr) to register any
>     YODA analysis object. But it's made a little fiddly because YODA returns
>     a stack-allocated object and what is needed is a new'd heap object. You
>     might need to do
> 
>     Scatter2D s = efficiency(*tmph1, *tmph2);
>     s.setPath(histoPath("foo"));
>     addAnalysisObject(s.newclone());
> 
>     or similar. It's not very pretty, which is why I should provide the
>     cosmetic wrapper functions.
> 
>     Andy
> 
>     --
>     Dr Andy Buckley, Royal Society University Research Fellow
>     Particle Physics Expt Group, University of Glasgow / PH Dept, CERN
> 
> 


-- 
Dr Andy Buckley, Royal Society University Research Fellow
Particle Physics Expt Group, University of Glasgow / PH Dept, CERN