Logo ROOT   6.10/09
Reference Guide
RooBifurGauss.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitModels *
4  * @(#)root/roofit:$Id$
5  * Authors: *
6  * Abi Soffer, Colorado State University, abi@slac.stanford.edu *
7  * *
8  * Copyright (c) 2000-2005, Regents of the University of California, *
9  * Colorado State University *
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 #include "RooFit.h"
17 
18 /** \class RooBifurGauss
19  \ingroup Roofit
20 
21 Bifurcated Gaussian p.d.f with different widths on left and right
22 side of maximum value
23 **/
24 
25 
26 #include "Riostream.h"
27 #include "TMath.h"
28 #include <math.h>
29 
30 #include "RooBifurGauss.h"
31 #include "RooAbsReal.h"
32 #include "RooMath.h"
33 
34 using namespace std;
35 
37 
38 ////////////////////////////////////////////////////////////////////////////////
39 
40 RooBifurGauss::RooBifurGauss(const char *name, const char *title,
41  RooAbsReal& _x, RooAbsReal& _mean,
42  RooAbsReal& _sigmaL, RooAbsReal& _sigmaR) :
43  RooAbsPdf(name, title),
44  x ("x" , "Dependent" , this, _x),
45  mean ("mean" , "Mean" , this, _mean),
46  sigmaL("sigmaL", "Left Sigma" , this, _sigmaL),
47  sigmaR("sigmaR", "Right Sigma", this, _sigmaR)
48 
49 {
50 }
51 
52 ////////////////////////////////////////////////////////////////////////////////
53 
54 RooBifurGauss::RooBifurGauss(const RooBifurGauss& other, const char* name) :
55  RooAbsPdf(other,name), x("x",this,other.x), mean("mean",this,other.mean),
56  sigmaL("sigmaL",this,other.sigmaL), sigmaR("sigmaR", this, other.sigmaR)
57 {
58 }
59 
60 ////////////////////////////////////////////////////////////////////////////////
61 
63  Double_t arg = x - mean;
64 
65  Double_t coef(0.0);
66 
67  if (arg < 0.0){
68  if (TMath::Abs(sigmaL) > 1e-30) {
69  coef = -0.5/(sigmaL*sigmaL);
70  }
71  } else {
72  if (TMath::Abs(sigmaR) > 1e-30) {
73  coef = -0.5/(sigmaR*sigmaR);
74  }
75  }
76 
77  return exp(coef*arg*arg);
78 }
79 
80 ////////////////////////////////////////////////////////////////////////////////
81 
82 Int_t RooBifurGauss::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
83 {
84  if (matchArgs(allVars,analVars,x)) return 1 ;
85  return 0 ;
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 
90 Double_t RooBifurGauss::analyticalIntegral(Int_t code, const char* rangeName) const
91 {
92  switch(code) {
93  case 1:
94  {
95  static Double_t root2 = sqrt(2.) ;
96  static Double_t rootPiBy2 = sqrt(atan2(0.0,-1.0)/2.0);
97 
98 // Double_t coefL(0.0), coefR(0.0);
99 // if (TMath::Abs(sigmaL) > 1e-30) {
100 // coefL = -0.5/(sigmaL*sigmaL);
101 // }
102 
103 // if (TMath::Abs(sigmaR) > 1e-30) {
104 // coefR = -0.5/(sigmaR*sigmaR);
105 // }
106 
107  Double_t xscaleL = root2*sigmaL;
108  Double_t xscaleR = root2*sigmaR;
109 
110  Double_t integral = 0.0;
111  if(x.max(rangeName) < mean)
112  {
113  integral = sigmaL * ( RooMath::erf((x.max(rangeName) - mean)/xscaleL) - RooMath::erf((x.min(rangeName) - mean)/xscaleL) );
114  }
115  else if (x.min(rangeName) > mean)
116  {
117  integral = sigmaR * ( RooMath::erf((x.max(rangeName) - mean)/xscaleR) - RooMath::erf((x.min(rangeName) - mean)/xscaleR) );
118  }
119  else
120  {
121  integral = sigmaR*RooMath::erf((x.max(rangeName) - mean)/xscaleR) - sigmaL*RooMath::erf((x.min(rangeName) - mean)/xscaleL);
122  }
123  // return rootPiBy2*(sigmaR*RooMath::erf((x.max(rangeName) - mean)/xscaleR) -
124  // sigmaL*RooMath::erf((x.min(rangeName) - mean)/xscaleL));
125  return integral*rootPiBy2;
126  }
127  }
128 
129  assert(0) ;
130  return 0 ; // to prevent compiler warnings
131 }
RooRealProxy x
Definition: RooBifurGauss.h:40
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
int Int_t
Definition: RtypesCore.h:41
RooRealProxy sigmaL
Definition: RooBifurGauss.h:42
STL namespace.
Short_t Abs(Short_t d)
Definition: TMathBase.h:108
Double_t evaluate() const
double sqrt(double)
Double_t x[n]
Definition: legend1.C:17
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
RooRealProxy sigmaR
Definition: RooBifurGauss.h:43
RooRealProxy mean
Definition: RooBifurGauss.h:41
Bifurcated Gaussian p.d.f with different widths on left and right side of maximum value...
Definition: RooBifurGauss.h:24
Double_t Mean(Long64_t n, const T *a, const Double_t *w=0)
Definition: TMath.h:973
#define ClassImp(name)
Definition: Rtypes.h:336
Double_t min(const char *rname=0) const
Definition: RooRealProxy.h:56
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)
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
static std::complex< double > erf(const std::complex< double > z)
complex erf function
Definition: RooMath.cxx:580
constexpr Double_t Sigma()
Definition: TMath.h:192
double exp(double)