Logo ROOT  
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 *
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 
17 /** \class RooGaussian
18  \ingroup Roofit
19 
20 Plain Gaussian p.d.f
21 **/
22 
23 #include "RooGaussian.h"
24 
25 #include "RooFit.h"
26 #include "BatchHelpers.h"
27 #include "RooRandom.h"
28 #include "RooMath.h"
29 #include "RooHelpers.h"
30 #include "RooVDTHeaders.h"
31 #include "RooFitComputeInterface.h"
32 
33 using namespace BatchHelpers;
34 using namespace std;
35 
37 
38 ////////////////////////////////////////////////////////////////////////////////
39 
40 RooGaussian::RooGaussian(const char *name, const char *title,
41  RooAbsReal& _x, RooAbsReal& _mean,
42  RooAbsReal& _sigma) :
43  RooAbsPdf(name,title),
44  x("x","Observable",this,_x),
45  mean("mean","Mean",this,_mean),
46  sigma("sigma","Width",this,_sigma)
47 {
48  RooHelpers::checkRangeOfParameters(this, {&_sigma}, 0);
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  const double arg = x - mean;
64  const double sig = sigma;
65  return exp(-0.5*arg*arg/(sig*sig));
66 }
67 
68 
69 ////////////////////////////////////////////////////////////////////////////////
70 /// Compute multiple values of Gaussian distribution.
72  return RooFitCompute::dispatch->computeGaussian(this, evalData, x->getValues(evalData, normSet), mean->getValues(evalData, normSet), sigma->getValues(evalData, normSet));
73 }
74 
75 ////////////////////////////////////////////////////////////////////////////////
76 
77 Int_t RooGaussian::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
78 {
79  if (matchArgs(allVars,analVars,x)) return 1 ;
80  if (matchArgs(allVars,analVars,mean)) return 2 ;
81  return 0 ;
82 }
83 
84 ////////////////////////////////////////////////////////////////////////////////
85 
86 Double_t RooGaussian::analyticalIntegral(Int_t code, const char* rangeName) const
87 {
88  assert(code==1 || code==2);
89 
90  //The normalisation constant 1./sqrt(2*pi*sigma^2) is left out in evaluate().
91  //Therefore, the integral is scaled up by that amount to make RooFit normalise
92  //correctly.
93  const double resultScale = sqrt(TMath::TwoPi()) * sigma;
94 
95  //Here everything is scaled and shifted into a standard normal distribution:
96  const double xscale = TMath::Sqrt2() * sigma;
97  double max = 0.;
98  double min = 0.;
99  if (code == 1){
100  max = (x.max(rangeName)-mean)/xscale;
101  min = (x.min(rangeName)-mean)/xscale;
102  } else { //No == 2 test because of assert
103  max = (mean.max(rangeName)-x)/xscale;
104  min = (mean.min(rangeName)-x)/xscale;
105  }
106 
107 
108  //Here we go for maximum precision: We compute all integrals in the UPPER
109  //tail of the Gaussian, because erfc has the highest precision there.
110  //Therefore, the different cases for range limits in the negative hemisphere are mapped onto
111  //the equivalent points in the upper hemisphere using erfc(-x) = 2. - erfc(x)
112  const double ecmin = std::erfc(std::abs(min));
113  const double ecmax = std::erfc(std::abs(max));
114 
115 
116  return resultScale * 0.5 * (
117  min*max < 0.0 ? 2.0 - (ecmin + ecmax)
118  : max <= 0. ? ecmax - ecmin : ecmin - ecmax
119  );
120 }
121 
122 ////////////////////////////////////////////////////////////////////////////////
123 
124 Int_t RooGaussian::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
125 {
126  if (matchArgs(directVars,generateVars,x)) return 1 ;
127  if (matchArgs(directVars,generateVars,mean)) return 2 ;
128  return 0 ;
129 }
130 
131 ////////////////////////////////////////////////////////////////////////////////
132 
134 {
135  assert(code==1 || code==2) ;
136  Double_t xgen ;
137  if(code==1){
138  while(1) {
140  if (xgen<x.max() && xgen>x.min()) {
141  x = xgen ;
142  break;
143  }
144  }
145  } else if(code==2){
146  while(1) {
148  if (xgen<mean.max() && xgen>mean.min()) {
149  mean = xgen ;
150  break;
151  }
152  }
153  } else {
154  cout << "error in RooGaussian generateEvent"<< endl;
155  }
156 
157  return;
158 }
RooGaussian::x
RooRealProxy x
Definition: RooGaussian.h:43
RooHelpers.h
TRandom::Gaus
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:274
RooFit.h
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
RooHelpers::checkRangeOfParameters
void checkRangeOfParameters(const RooAbsReal *callingClass, std::initializer_list< const RooAbsReal * > pars, double min=-std::numeric_limits< double >::max(), double max=std::numeric_limits< double >::max(), bool limitsInAllowedRange=false, std::string extraMessage="")
Check if the parameters have a range, and warn if the range extends below / above the set limits.
Definition: RooHelpers.cxx:117
RooGaussian::evaluate
Double_t evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
Definition: RooGaussian.cxx:61
RooGaussian.h
exp
double exp(double)
x
Double_t x[n]
Definition: legend1.C:17
RooGaussian
Plain Gaussian p.d.f.
Definition: RooGaussian.h:24
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooGaussian::analyticalIntegral
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const override
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
Definition: RooGaussian.cxx:86
RooFitComputeInterface.h
ROOT::Math::erfc
double erfc(double x)
Complementary error function.
Definition: SpecFuncMathCore.cxx:44
bool
BatchHelpers::RunContext
Data that has to be passed around when evaluating functions / PDFs.
Definition: RunContext.h:32
RooGaussian::getGenerator
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const override
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooGaussian.cxx:124
BatchHelpers
Definition: RooFitComputeInterface.h:9
RooGaussian::RooGaussian
RooGaussian()
Definition: RooGaussian.h:26
RooGaussian::evaluateSpan
RooSpan< double > evaluateSpan(BatchHelpers::RunContext &evalData, const RooArgSet *normSet) const override
Compute multiple values of Gaussian distribution.
Definition: RooGaussian.cxx:71
RooTemplateProxy::min
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
Definition: RooTemplateProxy.h:275
RooGaussian::generateEvent
void generateEvent(Int_t code) override
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooGaussian.cxx:133
RooRandom.h
RooVDTHeaders.h
RooAbsReal::getValues
virtual RooSpan< const double > getValues(BatchHelpers::RunContext &evalData, const RooArgSet *normSet=nullptr) const
by this change, please consult the release notes for ROOT 6.24 for guidance on how to make this trans...
Definition: RooAbsReal.cxx:313
RooFitCompute::RooFitComputeInterface::computeGaussian
virtual RooSpan< double > computeGaussian(const RooAbsReal *, BatchHelpers::RunContext &, RooSpan< const double > x, RooSpan< const double > mean, RooSpan< const double > sigma)=0
BatchHelpers.h
sqrt
double sqrt(double)
RooGaussian::sigma
RooRealProxy sigma
Definition: RooGaussian.h:45
TMath::TwoPi
constexpr Double_t TwoPi()
Definition: TMath.h:44
sigma
const Double_t sigma
Definition: h1analysisProxy.h:11
Double_t
double Double_t
Definition: RtypesCore.h:59
RooGaussian::mean
RooRealProxy mean
Definition: RooGaussian.h:44
RooFitCompute::dispatch
R__EXTERN RooFitComputeInterface * dispatch
The dispatch pointer always points to the most recently created RooFitComputeClass object.
Definition: RooFitComputeInterface.h:46
RooGaussian::getAnalyticalIntegral
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
Definition: RooGaussian.cxx:77
name
char name[80]
Definition: TGX11.cxx:110
TMath::Sqrt2
constexpr Double_t Sqrt2()
Definition: TMath.h:88
RooTemplateProxy::max
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
Definition: RooTemplateProxy.h:277
RooAbsReal::matchArgs
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Definition: RooAbsReal.cxx:3407
RooRandom::randomGenerator
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:53
RooAbsPdf
Definition: RooAbsPdf.h:43
RooMath.h
RooSpan
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
int