ROOT   Reference Guide
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 *
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
20RooArgusBG is a RooAbsPdf implementation describing the ARGUS background shape.
21\f[
22 \mathrm{Argus}(m, m_0, c, p) = \mathcal{N} \cdot m \cdot \left[ 1 - \left( \frac{m}{m_0} \right)^2 \right]^p
23 \cdot \exp\left[ c \cdot \left(1 - \left(\frac{m}{m_0}\right)^2 \right) \right]
24\f]
25\image html RooArgusBG.png
26*/
27
28#include "RooArgusBG.h"
29#include "RooRealVar.h"
30#include "RooRealConstant.h"
31#include "RooMath.h"
32#include "RooBatchCompute.h"
33
34#include "TMath.h"
35
36#include <cmath>
37using namespace std;
38
40
41////////////////////////////////////////////////////////////////////////////////
42/// Constructor.
43
44RooArgusBG::RooArgusBG(const char *name, const char *title,
45 RooAbsReal& _m, RooAbsReal& _m0, RooAbsReal& _c) :
46 RooAbsPdf(name, title),
47 m("m","Mass",this,_m),
48 m0("m0","Resonance mass",this,_m0),
49 c("c","Slope parameter",this,_c),
50 p("p","Power",this,(RooRealVar&)RooRealConstant::value(0.5))
51{
52}
53
54////////////////////////////////////////////////////////////////////////////////
55/// Constructor.
56
57RooArgusBG::RooArgusBG(const char *name, const char *title,
58 RooAbsReal& _m, RooAbsReal& _m0, RooAbsReal& _c, RooAbsReal& _p) :
59 RooAbsPdf(name, title),
60 m("m","Mass",this,_m),
61 m0("m0","Resonance mass",this,_m0),
62 c("c","Slope parameter",this,_c),
63 p("p","Power",this,_p)
64{
65}
66
67////////////////////////////////////////////////////////////////////////////////
68/// Constructor.
69
70RooArgusBG::RooArgusBG(const RooArgusBG& other, const char* name) :
71 RooAbsPdf(other,name),
72 m("m",this,other.m),
73 m0("m0",this,other.m0),
74 c("c",this,other.c),
75 p("p",this,other.p)
76{
77}
78
79////////////////////////////////////////////////////////////////////////////////
80
81double RooArgusBG::evaluate() const {
82 double t= m/m0;
83 if(t >= 1) return 0;
84
85 double u= 1 - t*t;
86 return m*TMath::Power(u,p)*exp(c*u) ;
87}
88
89////////////////////////////////////////////////////////////////////////////////
90
91void RooArgusBG::computeBatch(cudaStream_t* stream, double* output, size_t nEvents, RooFit::Detail::DataMap const& dataMap) const
92{
94 dispatch->compute(stream, RooBatchCompute::ArgusBG, output, nEvents,
95 {dataMap.at(m), dataMap.at(m0), dataMap.at(c), dataMap.at(p)});
96}
97
98////////////////////////////////////////////////////////////////////////////////
99
100Int_t RooArgusBG::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
101{
102 if (p.arg().isConstant()) {
103 // We can integrate over m if power = 0.5
104 if (matchArgs(allVars,analVars,m) && p == 0.5) return 1;
105 }
106 return 0;
107
108}
109
110////////////////////////////////////////////////////////////////////////////////
111
112double RooArgusBG::analyticalIntegral(Int_t code, const char* rangeName) const
113{
114 R__ASSERT(code==1);
115 // Formula for integration over m when p=0.5
116 static const double pi = atan2(0.0,-1.0);
117 double min = (m.min(rangeName) < m0) ? m.min(rangeName) : m0;
118 double max = (m.max(rangeName) < m0) ? m.max(rangeName) : m0;
119 double f1 = (1.-TMath::Power(min/m0,2));
120 double f2 = (1.-TMath::Power(max/m0,2));
121 double aLow, aHigh ;
122 if ( c < 0. ) {
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 } else if ( c == 0. ) {
126 aLow = -m0*m0/3.*f1*sqrt(f1);
127 aHigh = -m0*m0/3.*f1*sqrt(f2);
128 } else {
129 aLow = 0.5*m0*m0*exp(c*f1)/(c*sqrt(c)) * (0.5*sqrt(pi)*(RooMath::faddeeva(sqrt(c*f1))).imag() - sqrt(c*f1));
130 aHigh = 0.5*m0*m0*exp(c*f2)/(c*sqrt(c)) * (0.5*sqrt(pi)*(RooMath::faddeeva(sqrt(c*f2))).imag() - sqrt(c*f2));
131 }
132 double area = aHigh - aLow;
133 //cout << "c = " << c << "aHigh = " << aHigh << " aLow = " << aLow << " area = " << area << endl ;
134 return area;
135
136}
#define c(i)
Definition: RSha256.hxx:101
#define ClassImp(name)
Definition: Rtypes.h:375
#define R__ASSERT(e)
Definition: TError.h:118
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
char name[80]
Definition: TGX11.cxx:110
bool isConstant() const
Check if the "Constant" attribute is set.
Definition: RooAbsArg.h:377
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:56
RooArgusBG is a RooAbsPdf implementation describing the ARGUS background shape.
Definition: RooArgusBG.h:22
RooRealProxy m
Definition: RooArgusBG.h:37
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
Definition: RooArgusBG.cxx:81
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
Definition: RooArgusBG.cxx:112
RooRealProxy c
Definition: RooArgusBG.h:39
RooRealProxy p
Definition: RooArgusBG.h:40
void computeBatch(cudaStream_t *, double *output, size_t size, RooFit::Detail::DataMap const &) const override
Base function for computing multiple values of a RooAbsReal.
Definition: RooArgusBG.cxx:91
RooRealProxy m0
Definition: RooArgusBG.h:38
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
Definition: RooArgusBG.cxx:100
virtual void compute(cudaStream_t *, Computer, RestrictArr, size_t, const VarVector &, const ArgVector &={})=0
auto & at(RooAbsArg const *arg, RooAbsArg const *=nullptr)
Definition: DataMap.h:88
static std::complex< double > erf(const std::complex< double > z)
complex erf function
Definition: RooMath.cxx:60
static std::complex< double > faddeeva(std::complex< double > z)
evaluate Faddeeva function for complex argument
Definition: RooMath.cxx:31
RooRealConstant provides static functions to create and keep track of RooRealVar constants.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
double max(const char *rname=nullptr) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
const T & arg() const
Return reference to object held in proxy.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
RVec< PromoteTypes< T0, T1 > > atan2(const T0 &x, const RVec< T1 > &v)
Definition: RVec.hxx:1781
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1761
TF1 * f1
Definition: legend1.C:11
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
R__EXTERN RooBatchComputeInterface * dispatchCUDA
R__EXTERN RooBatchComputeInterface * dispatchCPU
This dispatch pointer points to an implementation of the compute library, provided one has been loade...
static constexpr double pi
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition: TMath.h:719
TMarker m
Definition: textangle.C:8
static void output()