[Rivet] CMS analysis of dijet ratios at large rapdities

[Vadim Oreshkin]

Dear Hendrik,

Here are the answers to your questions.

- In the main error formula "sqrt( ((1.-2.*w)*e1*e1 + w*w*e2*e2 )/(b2*b2))"
   the errors appear in fourth power, while in the ROOT version they
   don't.

The link to the source code is outdated (because the version of Root has changed).
For the current version of ROOT it should be
http://root.cern.ch/root/html/src/TH1.cxx.html#2746.

So the point is that the histogram division implemented in Rivet assumes
that the numerator and denominator
are statistically independent. However this is not the case for our analysis, because
the set of entries falling into the exclusive histogram is a subset of the set of entries
falling into the inclusive histogram. (If the number of exclusive entries changes,
then the number of inclusive entries changes too). The probability for an inclusive
entry to be an exclusive entry is equal to an inverse of our observable, which is called briefly "Kfactor".

That being said, the number of exclusive entries is distributed
according to binomial distribution with parameter p = 1/Kfactor.
(http://en.wikipedia.org/wiki/Binomial_distribution)
So we take the "width" of this binomial distribution as an error of number exclusive entries.
n(1-p)p, where n is an number of inclusive entries. After dividing by n*n (in order to obtain
the error for 1/Kfactor, rather than for n/Kfactor=number of exclusive entries), the n goes to denominator.

  fSumw2.fArray[bin] = TMath::Abs(w*(1-w)/b2);

For weighted versions of histograms we use a more complicated formula,
which we borrowed from ROOT:
http://root.cern.ch/phpBB3/viewtopic.php?f=3&t=3753


  fSumw2.fArray[bin] = TMath::Abs( ( (1-2*w)*e1*e1 + w*w*e2*e2 )/(b2*b2) );


- Numerator and denominator seem to be swapped compared to the ROOT
   version (b1 = h2.binHeight(i) and b2 = h1.binHeight(i)).
That's because our observable is a ratio of a set to a subset, rather than a subset to a set.
(This is indicated in the name of the function which makes the histogram division)
- At the end you divide the error by a factor of 2.

We divide by 2, because ROOT and Rivet use different conventions for storing errors.
Rivet setErrorPlus,setErrorMinus store the value given by the argument as
a pair assymmetric errors, while in ROOT we store single value of an error.

Best regards,
Vadim.


Dear All,

can one of you experts who did the analysis please reply to Hendrik and the rivet forum (rivet at project.hepforge.org<mailto:rivet at project.hepforge.org>). It is better that you discuss directly through the Rivet forum instead of using me as a middle hand.

Thanks and cheers,
Albert

ps. if you are not subscribed to the rivet forum you can do that here, if you want: http://www.hepforge.org/lists/listinfo/rivet




-------- Original Message --------
Subject:        Re: [Rivet] CMS analysis of dijet ratios at large rapdities
Date:   Mon, 1 Oct 2012 15:51:03 +0100
From:   Hendrik Hoeth <hendrik.hoeth at cern.ch><mailto:hendrik.hoeth at cern.ch>
To:     Albert Knutsson <albert.knutsson at desy.de><mailto:albert.knutsson at desy.de>
CC:     rivet at projects.hepforge.org<mailto:rivet at projects.hepforge.org>



Hi Albert,

thanks for the analysis. I've had a look at the code, and I have to
admit that I don't understand the error calculation in your
divide_set_by_subset_with_binomial_errors() function at all. There
is a reference to some ROOT code, but it points to an empty line
(http://root.cern.ch/root/html/src/TH1.cxx.html#2592) -- I can only
suspect that it's meant to point at the TH1::Divide(const TH1 *h1, const
TH1 *h2, Double_t c1, Double_t c2, Option_t *option) in line 2651 in the
current version of TH1.cxx.

Anyway, here is what I find suspicious:

- In the main error formula "sqrt( ((1.-2.*w)*e1*e1 + w*w*e2*e2 )/(b2*b2))"
  the errors appear in fourth power, while in the ROOT version they
  don't.

- Numerator and denominator seem to be swapped compared to the ROOT
  version (b1 = h2.binHeight(i) and b2 = h1.binHeight(i)).

- At the end you divide the error by a factor of 2.

Can you please confirm that this code is actually what you want to do?

Cheers,

    Hendrik

--
If your dreams don't scare you, then you are not dreaming big enough.