Logo ROOT   6.08/07
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 /**
17 \file RooBukinPdf.cxx
18 \class RooBukinPdf
19 \ingroup Roofit
20 
21 RooBukinPdf implements the NovosibirskA function
22 **/
23 
24 // Original Fortran Header below
25 /*****************************************************************************
26  * Fitting function for asymmetric peaks with 6 free parameters: *
27  * Ap - peak value *
28  * Xp - peak position *
29  * sigp - FWHM divided by 2*sqrt(2*log(2))=2.35 *
30  * xi - peak asymmetry parameter *
31  * rho1 - parameter of the "left tail" *
32  * rho2 - parameter of the "right tail" *
33  * --------------------------------------------- *
34  * May 26, 2003 *
35  * A.Bukin, Budker INP, Novosibirsk *
36  * Documentation: *
37  * http://www.slac.stanford.edu/BFROOT/www/Organization/CollabMtgs/2003/detJuly2003/Tues3a/bukin.ps
38  * ------------------------------------------- *
39  *****************************************************************************/
40 
41 #include "RooFit.h"
42 
43 #include <math.h>
44 
45 
46 #include "RooBukinPdf.h"
47 #include "RooRealVar.h"
48 #include "TMath.h"
49 
50 using namespace std;
51 
53 
54 
55 
56 ////////////////////////////////////////////////////////////////////////////////
57 
58 RooBukinPdf::RooBukinPdf(const char *name, const char *title,
59  RooAbsReal& _x, RooAbsReal& _Xp,
60  RooAbsReal& _sigp, RooAbsReal& _xi,
61  RooAbsReal& _rho1, RooAbsReal& _rho2) :
62  // The two addresses refer to our first dependent variable and
63  // parameter, respectively, as declared in the rdl file
64  RooAbsPdf(name, title),
65  x("x","x",this,_x),
66  Xp("Xp","Xp",this,_Xp),
67  sigp("sigp","sigp",this,_sigp),
68  xi("xi","xi",this,_xi),
69  rho1("rho1","rho1",this,_rho1),
70  rho2("rho2","rho2",this,_rho2)
71 {
72  // Constructor
73  consts = 2*sqrt(2*log(2.));
74 }
75 
76 
77 
78 ////////////////////////////////////////////////////////////////////////////////
79 
80 RooBukinPdf::RooBukinPdf(const RooBukinPdf& other, const char *name):
81  RooAbsPdf(other,name),
82  x("x",this,other.x),
83  Xp("Xp",this,other.Xp),
84  sigp("sigp",this,other.sigp),
85  xi("xi",this,other.xi),
86  rho1("rho1",this,other.rho1),
87  rho2("rho2",this,other.rho2)
88 
89 {
90  // Copy constructor
91  consts = 2*sqrt(2*log(2.));
92 }
93 
94 
95 
96 ////////////////////////////////////////////////////////////////////////////////
97 /// Implementation
98 
100 {
101  double r1=0,r2=0,r3=0,r4=0,r5=0,hp=0;
102  double x1 = 0,x2 = 0;
103  double fit_result = 0;
104 
105  hp=sigp*consts;
106  r3=log(2.);
107  r4=sqrt(TMath::Power(xi,2)+1);
108  r1=xi/r4;
109 
110  if(TMath::Abs(xi) > exp(-6.)){
111  r5=xi/log(r4+xi);
112  }
113  else
114  r5=1;
115 
116  x1 = Xp + (hp / 2) * (r1-1);
117  x2 = Xp + (hp / 2) * (r1+1);
118 
119  //--- Left Side
120  if(x < x1){
121  r2=rho1*TMath::Power((x-x1)/(Xp-x1),2)-r3 + 4 * r3 * (x-x1)/hp * r5 * r4/TMath::Power((r4-xi),2);
122  }
123 
124 
125  //--- Center
126  else if(x < x2) {
127  if(TMath::Abs(xi) > exp(-6.)) {
128  r2=log(1 + 4 * xi * r4 * (x-Xp)/hp)/log(1+2*xi*(xi-r4));
129  r2=-r3*(TMath::Power(r2,2));
130  }
131  else{
132  r2=-4*r3*TMath::Power(((x-Xp)/hp),2);
133  }
134  }
135 
136 
137  //--- Right Side
138  else {
139  r2=rho2*TMath::Power((x-x2)/(Xp-x2),2)-r3 - 4 * r3 * (x-x2)/hp * r5 * r4/TMath::Power((r4+xi),2);
140  }
141 
142  if(TMath::Abs(r2) > 100){
143  fit_result = 0;
144  }
145  else{
146  //---- Normalize the result
147  fit_result = exp(r2);
148  }
149 
150  return fit_result;
151 
152 }
RooRealProxy rho2
Definition: RooBukinPdf.h:65
RooRealProxy rho1
Definition: RooBukinPdf.h:64
RooRealProxy x
Definition: RooBukinPdf.h:60
STL namespace.
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:501
double sqrt(double)
static const double x2[5]
Double_t x[n]
Definition: legend1.C:17
unsigned int r3[N_CITIES]
Definition: simanTSP.cxx:323
RooRealProxy sigp
Definition: RooBukinPdf.h:62
RooRealProxy xi
Definition: RooBukinPdf.h:63
double consts
Definition: RooBukinPdf.h:71
unsigned int r1[N_CITIES]
Definition: simanTSP.cxx:321
static const double x1[5]
#define ClassImp(name)
Definition: Rtypes.h:279
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
RooRealProxy Xp
Definition: RooBukinPdf.h:61
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
Double_t evaluate() const
Implementation.
Definition: RooBukinPdf.cxx:99
RooBukinPdf implements the NovosibirskA function.
Definition: RooBukinPdf.h:45
double exp(double)
unsigned int r2[N_CITIES]
Definition: simanTSP.cxx:322
char name[80]
Definition: TGX11.cxx:109
double log(double)