Logo ROOT   6.18/05
Reference Guide
RooBukinPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * RW, Ruddick William UC Colorado wor@slac.stanford.edu *
7 * *
8 * Copyright (c) 2000-2005, Regents of the University of California *
9 * and Stanford University. All rights reserved. *
10 * *
11 * Redistribution and use in source and binary forms, *
12 * with or without modification, are permitted according to the terms *
13 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
14 *****************************************************************************/
15
16/** \class RooBukinPdf
17 \ingroup Roofit
18
19RooBukinPdf implements the NovosibirskA function. For the parameters, see
20RooBukinPdf().
21
22Credits:
23May 26, 2003.
24A.Bukin, Budker INP, Novosibirsk
25
26http://www.slac.stanford.edu/BFROOT/www/Organization/CollabMtgs/2003/detJuly2003/Tues3a/bukin.ps
27**/
28
29#include "RooFit.h"
30
31#include <math.h>
32
33
34#include "RooBukinPdf.h"
35#include "RooRealVar.h"
36#include "TMath.h"
37
38using namespace std;
39
41
42////////////////////////////////////////////////////////////////////////////////
43/// Construct a Bukin PDF.
44/// \param name The name of the PDF for RooFit's bookeeping.
45/// \param title The title for e.g. plotting it.
46/// \param _x The variable.
47/// \param _Xp The peak position.
48/// \param _sigp The peak width as FWHM divided by 2*sqrt(2*log(2))=2.35
49/// \param _xi Peak asymmetry. Use values around 0.
50/// \param _rho1 Left tail. Use slightly negative starting values.
51/// \param _rho2 Right tail. Use slightly positive starting values.
52RooBukinPdf::RooBukinPdf(const char *name, const char *title,
53 RooAbsReal& _x, RooAbsReal& _Xp,
54 RooAbsReal& _sigp, RooAbsReal& _xi,
55 RooAbsReal& _rho1, RooAbsReal& _rho2) :
56 // The two addresses refer to our first dependent variable and
57 // parameter, respectively, as declared in the rdl file
58 RooAbsPdf(name, title),
59 x("x","x",this,_x),
60 Xp("Xp","Xp",this,_Xp),
61 sigp("sigp","sigp",this,_sigp),
62 xi("xi","xi",this,_xi),
63 rho1("rho1","rho1",this,_rho1),
64 rho2("rho2","rho2",this,_rho2)
65{
66 // Constructor
67 consts = 2*sqrt(2*log(2.));
68}
69
70////////////////////////////////////////////////////////////////////////////////
71/// Copy a Bukin PDF.
72RooBukinPdf::RooBukinPdf(const RooBukinPdf& other, const char *name):
73 RooAbsPdf(other,name),
74 x("x",this,other.x),
75 Xp("Xp",this,other.Xp),
76 sigp("sigp",this,other.sigp),
77 xi("xi",this,other.xi),
78 rho1("rho1",this,other.rho1),
79 rho2("rho2",this,other.rho2)
80
81{
82 // Copy constructor
83 consts = 2*sqrt(2*log(2.));
84}
85
86////////////////////////////////////////////////////////////////////////////////
87/// Implementation
88
90{
91 double r1=0,r2=0,r3=0,r4=0,r5=0,hp=0;
92 double x1 = 0,x2 = 0;
93 double fit_result = 0;
94
95 hp=sigp*consts;
96 r3=log(2.);
97 r4=sqrt(TMath::Power(xi,2)+1);
98 r1=xi/r4;
99
100 if(TMath::Abs(xi) > exp(-6.)){
101 r5=xi/log(r4+xi);
102 }
103 else
104 r5=1;
105
106 x1 = Xp + (hp / 2) * (r1-1);
107 x2 = Xp + (hp / 2) * (r1+1);
108
109 //--- Left Side
110 if(x < x1){
111 r2=rho1*TMath::Power((x-x1)/(Xp-x1),2)-r3 + 4 * r3 * (x-x1)/hp * r5 * r4/TMath::Power((r4-xi),2);
112 }
113
114
115 //--- Center
116 else if(x < x2) {
117 if(TMath::Abs(xi) > exp(-6.)) {
118 r2=log(1 + 4 * xi * r4 * (x-Xp)/hp)/log(1+2*xi*(xi-r4));
119 r2=-r3*(TMath::Power(r2,2));
120 }
121 else{
122 r2=-4*r3*TMath::Power(((x-Xp)/hp),2);
123 }
124 }
125
126
127 //--- Right Side
128 else {
129 r2=rho2*TMath::Power((x-x2)/(Xp-x2),2)-r3 - 4 * r3 * (x-x2)/hp * r5 * r4/TMath::Power((r4+xi),2);
130 }
131
132 if(TMath::Abs(r2) > 100){
133 fit_result = 0;
134 }
135 else{
136 //---- Normalize the result
137 fit_result = exp(r2);
138 }
139
140 return fit_result;
141
142}
static const double x2[5]
static const double x1[5]
double Double_t
Definition: RtypesCore.h:55
#define ClassImp(name)
Definition: Rtypes.h:365
char name[80]
Definition: TGX11.cxx:109
double sqrt(double)
double exp(double)
double log(double)
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooBukinPdf implements the NovosibirskA function.
Definition: RooBukinPdf.h:29
RooRealProxy rho2
Definition: RooBukinPdf.h:49
RooRealProxy sigp
Definition: RooBukinPdf.h:46
double consts
Definition: RooBukinPdf.h:55
RooRealProxy rho1
Definition: RooBukinPdf.h:48
RooRealProxy x
Definition: RooBukinPdf.h:44
Double_t evaluate() const
Implementation.
Definition: RooBukinPdf.cxx:89
RooRealProxy xi
Definition: RooBukinPdf.h:47
RooRealProxy Xp
Definition: RooBukinPdf.h:45
Double_t x[n]
Definition: legend1.C:17
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:723
Short_t Abs(Short_t d)
Definition: TMathBase.h:120