ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
RooArgusBG.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitModels *
4  * @(#)root/roofit:$Id$
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 
17 /** \class RooArgusBG
18  \ingroup Roofit
19 
20 RooArgusBG is a RooAbsPdf implementation describing the ARGUS background shape.
21 */
22 
23 #include "RooFit.h"
24 
25 #include "Riostream.h"
26 #include "Riostream.h"
27 #include <math.h>
28 
29 #include "RooArgusBG.h"
30 #include "RooRealVar.h"
31 #include "RooRealConstant.h"
32 #include "RooMath.h"
33 #include "TMath.h"
34 
35 #include "TError.h"
36 
37 using namespace std;
38 
40 
41 
42 ////////////////////////////////////////////////////////////////////////////////
43 /// Constructor.
44 
45 RooArgusBG::RooArgusBG(const char *name, const char *title,
46  RooAbsReal& _m, RooAbsReal& _m0, RooAbsReal& _c) :
47  RooAbsPdf(name, title),
48  m("m","Mass",this,_m),
49  m0("m0","Resonance mass",this,_m0),
50  c("c","Slope parameter",this,_c),
51  p("p","Power",this,(RooRealVar&)RooRealConstant::value(0.5))
52 {
53 }
54 
55 
56 ////////////////////////////////////////////////////////////////////////////////
57 /// Constructor.
58 
59 RooArgusBG::RooArgusBG(const char *name, const char *title,
60  RooAbsReal& _m, RooAbsReal& _m0, RooAbsReal& _c, RooAbsReal& _p) :
61  RooAbsPdf(name, title),
62  m("m","Mass",this,_m),
63  m0("m0","Resonance mass",this,_m0),
64  c("c","Slope parameter",this,_c),
65  p("p","Power",this,_p)
66 {
67 }
68 
69 
70 ////////////////////////////////////////////////////////////////////////////////
71 /// Constructor.
72 
73 RooArgusBG::RooArgusBG(const RooArgusBG& other, const char* name) :
74  RooAbsPdf(other,name),
75  m("m",this,other.m),
76  m0("m0",this,other.m0),
77  c("c",this,other.c),
78  p("p",this,other.p)
79 {
80 }
81 
82 
83 
84 ////////////////////////////////////////////////////////////////////////////////
85 
87  Double_t t= m/m0;
88  if(t >= 1) return 0;
89 
90  Double_t u= 1 - t*t;
91  //cout << "c = " << c << " result = " << m*TMath::Power(u,p)*exp(c*u) << endl ;
92  return m*TMath::Power(u,p)*exp(c*u) ;
93 }
94 
95 
96 
97 ////////////////////////////////////////////////////////////////////////////////
98 
99 Int_t RooArgusBG::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
100 {
101  if (p.arg().isConstant()) {
102  // We can integrate over m if power = 0.5
103  if (matchArgs(allVars,analVars,m) && p == 0.5) return 1;
104  }
105  return 0;
106 
107 }
108 
109 
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 
113 Double_t RooArgusBG::analyticalIntegral(Int_t code, const char* rangeName) const
114 {
115  R__ASSERT(code==1);
116  // Formula for integration over m when p=0.5
117  static const Double_t pi = atan2(0.0,-1.0);
118  Double_t min = (m.min(rangeName) < m0) ? m.min(rangeName) : m0;
119  Double_t max = (m.max(rangeName) < m0) ? m.max(rangeName) : m0;
120  Double_t f1 = (1.-TMath::Power(min/m0,2));
121  Double_t f2 = (1.-TMath::Power(max/m0,2));
122  Double_t aLow, aHigh ;
123  aLow = -0.5*m0*m0*(exp(c*f1)*sqrt(f1)/c + 0.5/TMath::Power(-c,1.5)*sqrt(pi)*RooMath::erf(sqrt(-c*f1)));
124  aHigh = -0.5*m0*m0*(exp(c*f2)*sqrt(f2)/c + 0.5/TMath::Power(-c,1.5)*sqrt(pi)*RooMath::erf(sqrt(-c*f2)));
125  Double_t area = aHigh - aLow;
126  //cout << "c = " << c << "aHigh = " << aHigh << " aLow = " << aLow << " area = " << area << endl ;
127  return area;
128 
129 }
130 
131 
RooRealProxy m0
Definition: RooArgusBG.h:41
static Vc_ALWAYS_INLINE int_v min(const int_v &x, const int_v &y)
Definition: vector.h:433
const double pi
return c
tuple f2
Definition: surfaces.py:24
#define R__ASSERT(e)
Definition: TError.h:98
int Int_t
Definition: RtypesCore.h:41
RooRealConstant provides static functions to create and keep track of RooRealVar constants.
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:501
double sqrt(double)
RooArgusBG is a RooAbsPdf implementation describing the ARGUS background shape.
Definition: RooArgusBG.h:25
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
Definition: RooArgusBG.cxx:99
RooRealProxy c
Definition: RooArgusBG.h:42
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
TThread * t[5]
Definition: threadsh1.C:13
TPaveLabel title(3, 27.1, 15, 28.7,"ROOT Environment and Tools")
TMarker * m
Definition: textangle.C:8
ClassImp(RooArgusBG) RooArgusBG
Constructor.
Definition: RooArgusBG.cxx:39
Bool_t isConstant() const
Definition: RooAbsArg.h:266
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Definition: RooArgusBG.cxx:113
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
double atan2(double, double)
static Vc_ALWAYS_INLINE int_v max(const int_v &x, const int_v &y)
Definition: vector.h:440
#define name(a, b)
Definition: linkTestLib0.cpp:5
Double_t evaluate() const
Definition: RooArgusBG.cxx:86
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
Double_t min(const char *rname=0) const
Definition: RooRealProxy.h:56
RooRealProxy p
Definition: RooArgusBG.h:43
TF1 * f1
Definition: legend1.C:11
static std::complex< double > erf(const std::complex< double > z)
complex erf function
Definition: RooMath.cxx:584
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
double exp(double)
float value
Definition: math.cpp:443
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
RooRealProxy m
Definition: RooArgusBG.h:40