ROOT   6.10/09 Reference Guide
RooGaussian.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 *
15  *****************************************************************************/
16
17 /** \class RooGaussian
18  \ingroup Roofit
19
20 Plain Gaussian p.d.f
21 **/
22
23 #include "RooFit.h"
24
25 #include "Riostream.h"
26 #include "Riostream.h"
27 #include <math.h>
28
29 #include "RooGaussian.h"
30 #include "RooAbsReal.h"
31 #include "RooRealVar.h"
32 #include "RooRandom.h"
33 #include "RooMath.h"
34
35 using namespace std;
36
38
39 ////////////////////////////////////////////////////////////////////////////////
40
41 RooGaussian::RooGaussian(const char *name, const char *title,
42  RooAbsReal& _x, RooAbsReal& _mean,
43  RooAbsReal& _sigma) :
44  RooAbsPdf(name,title),
45  x("x","Observable",this,_x),
46  mean("mean","Mean",this,_mean),
47  sigma("sigma","Width",this,_sigma)
48 {
49 }
50
51 ////////////////////////////////////////////////////////////////////////////////
52
53 RooGaussian::RooGaussian(const RooGaussian& other, const char* name) :
54  RooAbsPdf(other,name), x("x",this,other.x), mean("mean",this,other.mean),
55  sigma("sigma",this,other.sigma)
56 {
57 }
58
59 ////////////////////////////////////////////////////////////////////////////////
60
62 {
63  Double_t arg= x - mean;
64  Double_t sig = sigma ;
65  Double_t ret =exp(-0.5*arg*arg/(sig*sig)) ;
66 // if (gDebug>2) {
67 // cout << "gauss(" << GetName() << ") x = " << x << " mean = " << mean << " sigma = " << sigma << " ret = " << ret << endl ;
68 // }
69  return ret ;
70 }
71
72 ////////////////////////////////////////////////////////////////////////////////
73 /// calculate and return the negative log-likelihood of the Poisson
74
76 {
77  return RooAbsPdf::getLogVal(set) ;
78 // Double_t prob = getVal(set) ;
79 // return log(prob) ;
80
81  Double_t arg= x - mean;
82  Double_t sig = sigma ;
83
84  //static const Double_t rootPiBy2 = sqrt(atan2(0.0,-1.0)/2.0);
85  //Double_t extra = -0.5*arg*arg/(sig*sig) - log(2*rootPiBy2*sig) ;
86  Double_t extra = -0.5*arg*arg/(sig*sig) - log(analyticalIntegral(1,0)) ;
87
88  return extra ;
89
90 }
91
92 ////////////////////////////////////////////////////////////////////////////////
93
94 Int_t RooGaussian::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
95 {
96  if (matchArgs(allVars,analVars,x)) return 1 ;
97  if (matchArgs(allVars,analVars,mean)) return 2 ;
98  return 0 ;
99 }
100
101 ////////////////////////////////////////////////////////////////////////////////
102
103 Double_t RooGaussian::analyticalIntegral(Int_t code, const char* rangeName) const
104 {
105  assert(code==1 || code==2) ;
106
107  static const Double_t root2 = sqrt(2.) ;
108  static const Double_t rootPiBy2 = sqrt(atan2(0.0,-1.0)/2.0);
109  Double_t xscale = root2*sigma;
110  Double_t ret = 0;
111  if(code==1){
112  ret = rootPiBy2*sigma*(RooMath::erf((x.max(rangeName)-mean)/xscale)-RooMath::erf((x.min(rangeName)-mean)/xscale));
113 // if (gDebug>2) {
114 // cout << "Int_gauss_dx(mean=" << mean << ",sigma=" << sigma << ", xmin=" << x.min(rangeName) << ", xmax=" << x.max(rangeName) << ")=" << ret << endl ;
115 // }
116  } else if(code==2) {
117  ret = rootPiBy2*sigma*(RooMath::erf((mean.max(rangeName)-x)/xscale)-RooMath::erf((mean.min(rangeName)-x)/xscale));
118  } else{
119  cout << "Error in RooGaussian::analyticalIntegral" << endl;
120  }
121  return ret ;
122
123 }
124
125 ////////////////////////////////////////////////////////////////////////////////
126
127 Int_t RooGaussian::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
128 {
129  if (matchArgs(directVars,generateVars,x)) return 1 ;
130  if (matchArgs(directVars,generateVars,mean)) return 2 ;
131  return 0 ;
132 }
133
134 ////////////////////////////////////////////////////////////////////////////////
135
137 {
138  assert(code==1 || code==2) ;
139  Double_t xgen ;
140  if(code==1){
141  while(1) {
143  if (xgen<x.max() && xgen>x.min()) {
144  x = xgen ;
145  break;
146  }
147  }
148  } else if(code==2){
149  while(1) {
151  if (xgen<mean.max() && xgen>mean.min()) {
152  mean = xgen ;
153  break;
154  }
155  }
156  } else {
157  cout << "error in RooGaussian generateEvent"<< endl;
158  }
159
160  return;
161 }
RooRealProxy x
Definition: RooGaussian.h:44
Double_t evaluate() const
Definition: RooGaussian.cxx:61
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition: TRandom.cxx:235
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
virtual Double_t getLogVal(const RooArgSet *set=0) const
Return the log of the current value with given normalization An error message is printed if the argum...
Definition: RooAbsPdf.cxx:606
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
STL namespace.
double sqrt(double)
Double_t x[n]
Definition: legend1.C:17
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:54
const Double_t sigma
Double_t getLogVal(const RooArgSet *set) const
calculate and return the negative log-likelihood of the Poisson
Definition: RooGaussian.cxx:75
RooRealProxy sigma
Definition: RooGaussian.h:46
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
RooRealProxy mean
Definition: RooGaussian.h:45
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)
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
Definition: RooGaussian.cxx:94
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
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
double exp(double)
void generateEvent(Int_t code)
Interface for generation of anan event using the algorithm corresponding to the specified code...
double log(double)