-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathRooArgusGenBG.cxx
105 lines (95 loc) · 3.84 KB
/
RooArgusGenBG.cxx
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
/*****************************************************************************
* Project: BaBar detector at the SLAC PEP-II B-factory
* Package: RooDKDalitz
* File: $Id: RooArgusGenBG.cc,v 1.2 2007/11/08 14:46:13 martinee Exp $
* Authors:
* WV, Wouter Verkerke, UC Santa Barbara, [email protected] *
* DK, David Kirkby, UC Irvine, [email protected] *
* FMV, Fernando Martinez-Vidal, IFIC-Valencia, [email protected] *
* History:
* - Apr 20: creation
* - Nov 11, 2007, adapted for D*->D0 bkg (lower threshold)
*****************************************************************************/
//
// Compile and load the class using:
// root> .L RooArgusGenBG.cxx+
//
#include "RooFit.h"
#include "Riostream.h"
#include <math.h>
#include "RooMath.h"
#include "RooRealConstant.h"
#include "RooRealVar.h"
#include "TMath.h"
#include "TRegexp.h"
#include "RooArgusGenBG.h"
ClassImp(RooArgusGenBG)
RooArgusGenBG::RooArgusGenBG(const char *name,
const char *title,
RooAbsReal &_m,
RooAbsReal &_m0,
RooAbsReal &_c,
ThresholdType _thresholdType)
: RooAbsPdf(name, title)
, m("m", "Mass", this, _m)
, m0("m0", "Resonance mass", this, _m0)
, c("c", "Slope parameter", this, _c)
, c2("c2", "Slope parameter2", this, (RooRealVar &)RooRealConstant::value(0.))
, p("p", "Power", this, (RooRealVar &)RooRealConstant::value(0.5))
, thresholdType(_thresholdType) {}
RooArgusGenBG::RooArgusGenBG(const char *name,
const char *title,
RooAbsReal &_m,
RooAbsReal &_m0,
RooAbsReal &_c,
RooAbsReal &_p,
ThresholdType _thresholdType)
: RooAbsPdf(name, title)
, m("m", "Mass", this, _m)
, m0("m0", "Resonance mass", this, _m0)
, c("c", "Slope parameter", this, _c)
, c2("c2", "Slope parameter2", this, (RooRealVar &)RooRealConstant::value(0.))
, p("p", "Power", this, _p)
, thresholdType(_thresholdType) {}
RooArgusGenBG::RooArgusGenBG(const char *name,
const char *title,
RooAbsReal &_m,
RooAbsReal &_m0,
RooAbsReal &_c,
RooAbsReal &_c2,
RooAbsReal &_p,
ThresholdType _thresholdType)
: RooAbsPdf(name, title)
, m("m", "Mass", this, _m)
, m0("m0", "Resonance mass", this, _m0)
, c("c", "Slope parameter", this, _c)
, c2("c2", "Slope parameter2", this, _c2)
, p("p", "Power", this, _p)
, thresholdType(_thresholdType) {}
RooArgusGenBG::RooArgusGenBG(const RooArgusGenBG &other, const char *name)
: RooAbsPdf(other, name)
, m("m", this, other.m)
, m0("m0", this, other.m0)
, c("c", this, other.c)
, c2("c2", this, other.c2)
, p("p", this, other.p)
, thresholdType(other.thresholdType) {}
Double_t RooArgusGenBG::evaluate() const {
if(thresholdType != BrianLowerThreshold) {
Double_t t = m / m0;
if((t >= 1 && thresholdType == UpperThreshold) || (t <= 1 && thresholdType == LowerThreshold))
return 0;
Double_t u = thresholdType == UpperThreshold ? 1 - t * t : t * t - 1;
// cout << "c = " << c << " result = " << m*TMath::Power(u,p)*exp(c*u) << endl ;
assert(u >= 0);
return m * TMath::Power(u, p) * exp(c * u);
} else {
// u^a e^{-b*x+c*x^2)
Double_t t = m / m0;
if(t <= 1)
return 0;
Double_t u = t * t - 1;
assert(u >= 0);
return m * TMath::Power(u, p) * exp(c * u + c2 * u * u);
}
}