Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
RooHistError.h
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitCore *
4  * File: $Id: RooHistError.h,v 1.14 2007/05/11 09:11:30 verkerke Exp $
5  * Authors: *
6  * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7  * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8  * *
9  * Copyright (c) 2000-2005, Regents of the University of California *
10  * and Stanford University. All rights reserved. *
11  * *
12  * Redistribution and use in source and binary forms, *
13  * with or without modification, are permitted according to the terms *
14  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15  *****************************************************************************/
16 #ifndef ROO_HIST_ERROR
17 #define ROO_HIST_ERROR
18 
19 #include "Rtypes.h"
20 #include "RooNumber.h"
21 #include "RooAbsFunc.h"
22 #include <cmath>
23 #include <iostream>
24 
25 class RooHistError {
26 public:
27  static const RooHistError &instance();
28  virtual ~RooHistError() {} ;
29 
30  Bool_t getPoissonInterval(Int_t n, Double_t &mu1, Double_t &mu2, Double_t nSigma= 1) const;
31  Bool_t getBinomialIntervalAsym(Int_t n, Int_t m, Double_t &a1, Double_t &a2, Double_t nSigma= 1) const;
32  Bool_t getBinomialIntervalEff(Int_t n, Int_t m, Double_t &a1, Double_t &a2, Double_t nSigma= 1) const;
33  Bool_t getInterval(const RooAbsFunc *Qu, const RooAbsFunc *Ql, Double_t pointEstimate, Double_t stepSize,
34  Double_t &lo, Double_t &hi, Double_t nSigma) const;
35 
36  static RooAbsFunc *createPoissonSum(Int_t n) ;
37  static RooAbsFunc *createBinomialSum(Int_t n, Int_t m, Bool_t eff) ;
38 
39 private:
40 
41 
42  Bool_t getPoissonIntervalCalc(Int_t n, Double_t &mu1, Double_t &mu2, Double_t nSigma= 1) const;
43  Double_t _poissonLoLUT[1000] ;
44  Double_t _poissonHiLUT[1000] ;
45 
46  RooHistError();
47  Double_t seek(const RooAbsFunc &f, Double_t startAt, Double_t step, Double_t value) const;
48 
49  // -----------------------------------------------------------
50  // Define a 1-dim RooAbsFunc of mu that evaluates the sum:
51  //
52  // Q(n|mu) = Sum_{k=0}^{n} P(k|mu)
53  //
54  // where P(n|mu) = exp(-mu) mu**n / n! is the Poisson PDF.
55  // -----------------------------------------------------------
56  class PoissonSum : public RooAbsFunc {
57  public:
58  inline PoissonSum(Int_t n) : RooAbsFunc(1), _n(n) { }
59  inline Double_t operator()(const Double_t xvec[]) const {
60  Double_t mu(xvec[0]),result(1),factorial(1);
61  for(Int_t k= 1; k <= _n; k++) {
62  factorial*= k;
63  result+= pow(mu,k)/factorial;
64  }
65  return exp(-mu)*result;
66  };
67  inline Double_t getMinLimit(UInt_t /*index*/) const { return 0; }
68  inline Double_t getMaxLimit(UInt_t /*index*/) const { return RooNumber::infinity() ; }
69  private:
70  Int_t _n;
71  };
72 
73  // -----------------------------------------------------------
74  // Define a 1-dim RooAbsFunc of a that evaluates the sum:
75  //
76  // Q(n|n+m,a) = Sum_{k=0}^{n} B(k|n+m,a)
77  //
78  // where B(n|n+m,a) = (n+m)!/(n!m!) ((1+a)/2)**n ((1-a)/2)**m
79  // is the Binomial PDF.
80  // -----------------------------------------------------------
81  class BinomialSumAsym : public RooAbsFunc {
82  public:
83  BinomialSumAsym(Int_t n, Int_t m) : RooAbsFunc(1), _n1(n), _N1(n+m) {
84  }
85  inline Double_t operator()(const Double_t xvec[]) const
86  {
87  Double_t p1(0.5*(1+xvec[0])),p2(1-p1),result(0),fact1(1),fact2(1);
88  for(Int_t k= 0; k <= _n1; k++) {
89  if(k > 0) { fact2*= k; fact1*= _N1-k+1; }
90  result+= fact1/fact2*pow(p1,k)*pow(p2,_N1-k);
91  }
92  return result;
93  };
94 
95  inline Double_t getMinLimit(UInt_t /*index*/) const { return -1; }
96  inline Double_t getMaxLimit(UInt_t /*index*/) const { return +1; }
97 
98  private:
99  Int_t _n1 ; // WVE Solaris CC5 doesn't want _n or _N here (likely compiler bug)
100  Int_t _N1 ;
101  } ;
102 
103 
104  // -----------------------------------------------------------
105  // Define a 1-dim RooAbsFunc of a that evaluates the sum:
106  //
107  // Q(n|n+m,a) = Sum_{k=0}^{n} B(k|n+m,a)
108  //
109  // where B(n|n+m,a) = (n+m)!/(n!m!) ((1+a)/2)**n ((1-a)/2)**m
110  // is the Binomial PDF.
111  // -----------------------------------------------------------
112  class BinomialSumEff : public RooAbsFunc {
113  public:
114  BinomialSumEff(Int_t n, Int_t m) : RooAbsFunc(1), _n1(n), _N1(n+m) {
115  }
116  inline Double_t operator()(const Double_t xvec[]) const
117  {
118  Double_t p1(xvec[0]),p2(1-p1),result(0),fact1(1),fact2(1);
119  for(Int_t k= 0; k <= _n1; k++) {
120  if(k > 0) { fact2*= k; fact1*= _N1-k+1; }
121  result+= fact1/fact2*pow(p1,k)*pow(p2,_N1-k);
122  }
123  return result;
124  };
125 
126  inline Double_t getMinLimit(UInt_t /*index*/) const { return 0; }
127  inline Double_t getMaxLimit(UInt_t /*index*/) const { return +1; }
128 
129  private:
130  Int_t _n1 ; // WVE Solaris CC5 doesn't want _n or _N here (likely compiler bug)
131  Int_t _N1 ;
132  } ;
133 
134  ClassDef(RooHistError,1) // Utility class for calculating histogram errors
135 };
136 
137 #endif