51 mean(
"mean",
"Mean",this,_mean),
52 sigma(
"sigma",
"Width",this,_sigma),
66 mean(
"mean",
"Mean",this,_mean),
67 sigma(
"sigma",
"Width",this,_sigma),
68 msf(
"msf",
"Mean Scale Factor",this,_msSF),
69 ssf(
"ssf",
"Sigma Scale Factor",this,_msSF)
81 mean(
"mean",
"Mean",this,_mean),
82 sigma(
"sigma",
"Width",this,_sigma),
83 msf(
"msf",
"Mean Scale Factor",this,_meanSF),
84 ssf(
"ssf",
"Sigma Scale Factor",this,_sigmaSF)
92 _flatSFInt(other._flatSFInt),
93 _asympInt(other._asympInt),
94 mean(
"mean",this,other.mean),
96 msf(
"msf",this,other.msf),
97 ssf(
"ssf",this,other.ssf)
138 double param1 = arg1 ? arg1->
getVal() : 0.0;
139 double param2 = arg2 ? arg2->getVal() : 0.0;
146 auto xVals = dataMap.
at(
x);
147 auto meanVals = dataMap.
at(
mean);
148 auto meanSfVals = dataMap.
at(
msf);
149 auto sigmaVals = dataMap.
at(
sigma);
150 auto sigmaSfVals = dataMap.
at(
ssf);
154 const double zeroVal = 0.0;
159 double basisSign =
_basisCode - 10 * (basisType - 1) - 2;
170 {xVals, meanVals, meanSfVals, sigmaVals, sigmaSfVals, param1Vals}, extraArgs);
175 if (xVals.size() !=
size || meanVals.size() != 1 || meanSfVals.size() != 1 || sigmaVals.size() != 1 ||
176 sigmaSfVals.size() != 1 || param1Vals.size() != 1 || param2Vals.size() != 1) {
180 for (
unsigned int i = 0; i <
size; ++i) {
181 output[i] =
evaluate(xVals[i], meanVals[0] * meanSfVals[0], sigmaVals[0] * sigmaSfVals[0], param1Vals[0],
189 static double root2(std::sqrt(2.)) ;
190 static double root2pi(std::sqrt(2.*std::atan2(0.,-1.))) ;
191 static double rootpi(std::sqrt(std::atan2(0.,-1.))) ;
198 double dGamma = param2;
199 if (dGamma==0) basisType =
expBasis;
204 double result = std::exp(-0.5*xprime*xprime)/(
sigma*root2pi) ;
217 double _x = omega *tau ;
218 double _y = tau*dgamma/2;
219 double xprime = (
x-
mean)/tau ;
220 double c =
sigma/(root2*tau) ;
221 double u = xprime/(2*
c) ;
225 if (basisSign!=
Minus)
result += evalCerf(0,-u,
c).real();
226 if (basisSign!=
Plus)
result += evalCerf(0, u,
c).real();
233 if (_x==0.)
return result ;
234 if (basisSign!=
Minus)
result += -evalCerf(-_x,-u,
c).imag();
235 if (basisSign!=
Plus)
result += -evalCerf( _x, u,
c).imag();
242 if (basisSign!=
Minus)
result += evalCerf(-_x,-u,
c).real();
243 if (basisSign!=
Plus)
result += evalCerf( _x, u,
c).real();
250 int sgn = ( basisType ==
coshBasis ? +1 : -1 );
251 if (basisSign!=
Minus)
result += 0.5*( evalCerf(0,-u,
c*(1-_y)).real()+sgn*evalCerf(0,-u,
c*(1+_y)).real()) ;
252 if (basisSign!=
Plus)
result += 0.5*(sgn*evalCerf(0, u,
c*(1-_y)).real()+ evalCerf(0, u,
c*(1+_y)).real()) ;
261 double f1 = std::exp(-u*u);
262 return (xprime - 2*
c*
c)*f0 + (2*
c/rootpi)*
f1 ;
270 double f1 = std::exp(-u*u);
271 double x2c2 = xprime - 2*
c*
c;
272 return ( x2c2*x2c2*f0 + (2*
c/rootpi)*x2c2*
f1 + 2*
c*
c*f0 );
334 static const double root2 = std::sqrt(2.) ;
336 static const double rootpi = std::sqrt(std::atan2(0.0,-1.0));
341 if (code==2) ssfInt = (
ssf.
max(rangeName)-
ssf.
min(rangeName)) ;
350 if (dGamma==0) basisType =
expBasis;
354 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 1st form" << endl ;
356 double xpmin = (
x.
min(rangeName)-(
mean*
msf))/xscale ;
357 double xpmax = (
x.
max(rangeName)-(
mean*
msf))/xscale ;
378 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 2nd form" << endl ;
386 double umin = xpmin/(2*
c) ;
387 double umax = xpmax/(2*
c) ;
390 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 3d form tau=" << tau << endl ;
399 double _x = omega * tau ;
400 double _y = tau*dgamma/2;
403 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 4th form omega = " << omega <<
", tau = " << tau << endl ;
405 if (_x==0)
return result*ssfInt ;
414 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 5th form omega = " << omega <<
", tau = " << tau << endl ;
425 if (
verboseEval()>0) {cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 8th form tau=" << tau << endl ; }
427 int sgn = ( basisType ==
coshBasis ? +1 : -1 );
428 if (basisSign!=
Minus)
result += 0.5*(
evalCerfInt(+1,0,tau/(1-_y),-umin,-umax,
c*(1-_y)).real()+ sgn*
evalCerfInt(+1,0,tau/(1+_y),-umin,-umax,
c*(1+_y)).real());
429 if (basisSign!=
Plus)
result += 0.5*(sgn*
evalCerfInt(-1,0,tau/(1-_y), umin, umax,
c*(1-_y)).real()+
evalCerfInt(-1,0,tau/(1+_y), umin, umax,
c*(1+_y)).real());
436 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 6th form tau=" << tau << endl ;
439 double f1 = std::exp(-umax*umax) - std::exp(-umin*umin);
444 double f2 = tmp1 - tmp2;
445 double f3 = xpmax*tmp1 - xpmin*tmp2;
447 double expc2 = std::exp(
c*
c);
451 (1 - 2*
c*
c)*expc2*f2 +
458 if (
verboseEval()>0) cout <<
"RooGaussModel::analyticalIntegral(" <<
GetName() <<
") 7th form tau=" << tau << endl ;
462 double tmpA1 = std::exp(-umax*umax);
463 double tmpA2 = std::exp(-umin*umin);
465 double f1 = tmpA1 - tmpA2;
466 double f2 = umax*tmpA1 - umin*tmpA2;
471 double f3 = tmpB1 - tmpB2;
472 double f4 = xpmax*tmpB1 - xpmin*tmpB2;
473 double f5 = xpmax*xpmax*tmpB1 - xpmin*xpmin*tmpB2;
475 double expc2 = std::exp(
c*
c);
478 (4*
c/rootpi)*((1-
c*
c)*
f1 +
c*f2) +
479 (2*
c*
c*(2*
c*
c-1) + 2)*expc2*f3 - (4*
c*
c-2)*expc2*f4 + expc2*f5
491 std::complex<double> diff(2., 0.);
493 diff = evalCerf(_x,umin,
c);
494 diff -= evalCerf(_x,umax,
c);
498 diff *= std::complex<double>(1., _x);
499 diff *= tau / (1.+_x*_x);
507 if (
matchArgs(directVars,generateVars,
x))
return 1 ;
521 if (xgen<xmax && xgen>
xmin) {
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
virtual void computeBatch(cudaStream_t *, double *output, size_t size, RooFit::Detail::DataMap const &) const
Base function for computing multiple values of a RooAbsReal.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
RooSpan< const double > at(RooAbsArg const *arg, RooAbsArg const *caller=nullptr)
Class RooGaussModel implements a RooResolutionModel that models a Gaussian distribution.
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
std::complex< double > evalCerfInt(double sign, double wt, double tau, double umin, double umax, double c) const
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
void computeBatch(cudaStream_t *, double *output, size_t size, RooFit::Detail::DataMap const &) const override
Base function for computing multiple values of a RooAbsReal.
~RooGaussModel() override
Destructor.
void generateEvent(Int_t code) override
Interface for generation of an event using the algorithm corresponding to the specified code.
Int_t basisCode(const char *name) const override
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const override
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
double analyticalIntegral(Int_t code, const char *rangeName) const override
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
static BasisType getBasisType(int basisCode)
static std::complex< double > erfc(const std::complex< double > z)
complex erfc function
static std::complex< double > erf(const std::complex< double > z)
complex erf function
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
RooRealConstant provides static functions to create and keep track of RooRealVar constants.
RooRealVar represents a variable that can be changed from the outside.
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
Int_t _basisCode
Identifier code for selected basis function.
RooAbsRealLValue & convVar() const
Return the convolution variable of the resolution model.
const RooFormulaVar & basis() const
RooTemplateProxy< RooAbsRealLValue > x
Dependent/convolution variable.
A simple container to hold a batch of data values.
double max(const char *rname=nullptr) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
const T & arg() const
Return reference to object held in proxy.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
const char * GetName() const override
Returns name of object.
This is the base class for the ROOT Random number generators.
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...
R__EXTERN RooBatchComputeInterface * dispatchCUDA
R__EXTERN RooBatchComputeInterface * dispatchCPU
This dispatch pointer points to an implementation of the compute library, provided one has been loade...
std::vector< double > ArgVector
__roohost__ __roodevice__ STD::complex< double > evalCerfApprox(double _x, double u, double c)
use the approximation: erf(z) = exp(-z*z)/(STD::sqrt(pi)*z) to explicitly cancel the divergent exp(y*...
__roohost__ __roodevice__ STD::complex< double > evalCerf(double swt, double u, double c)