[Rivet] Efficiency histograms in Rivet 2

[James Robinson]

Dear Andy and David,

I've done a bit more digging here and the problematic check is on the
effective number of entries.
​I have two Histo1DPtr objects,
hist_denominator
and
hist_numerator
and ​
I'm filling both histograms with

hist
​->​
fill( value, weight );

where sometimes I fill hist_denominator and sometimes hist_numerator AND
hist_denominator
​. If I print out every time a particular bin gets filled, you can see
this.​


Rivet.Analysis.ATLAS_STDM_2012_17: INFO  Filling hist_denominator with
​weight ​
587738 at value = 3.00608
Rivet.Analysis.ATLAS_STDM_2012_17: INFO    bin now contains actual entries:
4 / 16
Rivet.Analysis.ATLAS_STDM_2012_17: INFO    bin now contains effective
entries: 1.15024 / 1.195 => 0.962543
Rivet.Analysis.ATLAS_STDM_2012_17: INFO  Filling hist_denominator
with weight 107673 at
​value
 = 3.40704
Rivet.Analysis.ATLAS_STDM_2012_17: INFO  Filling hist_numerator with weight
​ ​
107673 at value = 3.40704
Rivet.Analysis.ATLAS_STDM_2012_17: INFO    bin now contains actual entries:
5 / 17
Rivet.Analysis.ATLAS_STDM_2012_17: INFO    bin now contains effective
entries: 1.20029 / 1.19789 => 1.002

​The second fill here is the problematic one, causing the efficiency to go
above 1.0. Adding a new value to both the numerator and denominator
increases the number of effective entries in both cases (as would be
expected) but increases effNumEntries in the numerator by more than the
denominator. This is what's causing the problem.

​I don't think the logic in using sum(w)^2 / sum(w^2)​ to calculate
nEntries is wrong, but I'm not sure exactly what's going on here.
​
Any ideas
​?

I can give a more detailed breakdown of exactly what I fill these
histograms with if needed (as you can see it's only 17 events).

Thanks,
James​

On 23 September 2014 20:16, James Robinson <james.robinson at cern.ch> 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> 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>> 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>>>
>> 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
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://www.hepforge.org/lists-archive/rivet/attachments/20140930/9db79b73/attachment.html>