53 sigma(
"sigma",
"Width",this,_sigma),
54 rlife(
"rlife",
"Life time",this,_rlife),
57 _flip(type==Flipped),_nlo(nlo), _flatSFInt(
kFALSE), _asympInt(
kFALSE)
70 sigma(
"sigma",
"Width",this,_sigma),
71 rlife(
"rlife",
"Life time",this,_rlife),
72 ssf(
"ssf",
"Sigma Scale Factor",this,_rsSF),
73 rsf(
"rsf",
"RLife Scale Factor",this,_rsSF),
90 sigma(
"sigma",
"Width",this,_sigma),
91 rlife(
"rlife",
"Life time",this,_rlife),
92 ssf(
"ssf",
"Sigma Scale Factor",this,_sigmaSF),
93 rsf(
"rsf",
"RLife Scale Factor",this,_rlifeSF),
107 sigma(
"sigma",this,other.sigma),
108 rlife(
"rlife",this,other.rlife),
109 ssf(
"ssf",this,other.ssf),
110 rsf(
"rsf",this,other.rsf),
113 _flatSFInt(other._flatSFInt),
114 _asympInt(other._asympInt)
173 if (dGamma==0) basisType =
expBasis;
178 if (
verboseEval()>2) cout <<
"RooGExpModel::evaluate(" <<
GetName() <<
") 1st form" << endl ;
180 Double_t expArg = sig*sig/(2*rtau*rtau) + fsign*
x/rtau ;
184 result = 1/(2*rtau) *
exp(expArg) *
RooMath::erfc(sig/(root2*rtau) + fsign*
x/(root2*sig));
188 result = 1/(2*rtau) *
exp(expArg +
logErfC(sig/(root2*rtau) + fsign*
x/(root2*sig))) ;
208 if (
verboseEval()>2) cout <<
"RooGExpModel::evaluate(" <<
GetName() <<
") 2nd form" << endl ;
216 if (
verboseEval()>2) cout <<
"RooGExpModel::evaluate(" <<
GetName() <<
") 3d form tau=" << tau << endl ;
227 if (
verboseEval()>2) cout <<
"RooGExpModel::evaluate(" <<
GetName() <<
") 4th form omega = "
228 << omega <<
", tau = " << tau << endl ;
230 if (wt==0.)
return result ;
231 if (basisSign!=
Minus) result += -1*
calcSinConv(+1,sig,tau,omega,rtau,fsign).imag() ;
232 if (basisSign!=
Plus) result += -1*
calcSinConv(-1,sig,tau,omega,rtau,fsign).imag() ;
240 <<
") 5th form omega = " << omega <<
", tau = " << tau << endl ;
242 if (basisSign!=
Minus) result +=
calcSinConv(+1,sig,tau,omega,rtau,fsign).real() ;
243 if (basisSign!=
Plus) result +=
calcSinConv(-1,sig,tau,omega,rtau,fsign).real() ;
254 <<
") 6th form = " << dgamma <<
", tau = " << tau << endl;
259 Double_t tau1 = 1/(1/tau-dgamma/2) ;
260 Double_t tau2 = 1/(1/tau+dgamma/2) ;
274 <<
") 7th form = " << dgamma <<
", tau = " << tau << endl;
279 Double_t tau1 = 1/(1/tau-dgamma/2) ;
280 Double_t tau2 = 1/(1/tau+dgamma/2) ;
304 ans=
log(t)+(-z*z-1.26551223+t*(1.00002368+t*(0.37409196+t*(0.09678418+t*(-0.18628806+
305 t*(0.27886807+t*(-1.13520398+t*(1.48851587+t*(-0.82215223+t*0.17087277)))))))));
307 ans=
log(2.0-t*
exp(-z*z-1.26551223+t*(1.00002368+t*(0.37409196+t*(0.09678418+t*(-0.18628806+
308 t*(0.27886807+t*(-1.13520398+t*(1.48851587+t*(-0.82215223+t*0.17087277))))))))));
330 std::complex<Double_t> eins(1,0);
331 std::complex<Double_t> k(1/tau,sign*omega);
334 return (
evalCerf(-sign*omega*tau,u1,c1)+std::complex<Double_t>(
evalCerfRe(u2,c2),0)*fsign*sign) / (eins + k*fsign*sign*rtau) ;
385 if ((
sign<0)&&(
fabs(tau-rtau)<tau/260)) {
388 if (
fabs(xp/MeanTau)>300) {
392 cFly=1./(MeanTau*MeanTau*root2pi) *
393 exp(-(-xp/MeanTau-sig*sig/(2*MeanTau*MeanTau)))
395 -(sig*sig/MeanTau+xp)*(rootpi/root2)*
RooMath::erfc(sig/(root2*MeanTau)+xp/(root2*sig)));
399 Double_t a=sig/(root2*MeanTau)+xp/(root2*sig);
400 cFly += 1./(MeanTau*MeanTau)
401 *
exp(-(-xp/MeanTau-sig*sig/(2*MeanTau*MeanTau)))
402 *0.5/MeanTau*epsilon*epsilon*
403 (
exp(-a*a)*(sig/MeanTau*root2/rootpi
404 -(4*a*sig*sig)/(2*rootpi*MeanTau*MeanTau)
405 +(-4/rootpi+8*a*a/rootpi)/6
407 +2/rootpi*(sig*sig/(MeanTau*MeanTau)+xp/MeanTau)*
408 (sig/(root2*MeanTau)-a*(sig*sig)/(2*MeanTau*MeanTau))
409 +2/rootpi*((3*sig*sig)/(2*MeanTau*MeanTau)+xp/MeanTau+
410 0.5*
TMath::Power(sig*sig/(MeanTau*MeanTau)+xp/MeanTau,2))*sig/(root2*MeanTau))
411 -(2*sig*sig/(MeanTau*MeanTau)+xp/MeanTau+(sig*sig/(MeanTau*MeanTau)+xp/MeanTau)*
412 (3*sig*sig/(2*MeanTau*MeanTau)+xp/MeanTau)
419 Double_t expArg2 = sig*sig/(2*rtau*rtau)+xp/rtau ;
425 term1 =
exp(expArg1+logErfC(sig/(root2*tau)-
sign*xp/(root2*sig))) ; ;
430 term2 =
exp(expArg2+logErfC(sig/(root2*rtau)+xp/(root2*sig))) ; ;
433 cFly=(term1+
sign*term2)/(2*(tau+
sign*rtau));
451 Double_t RooGExpModel::calcCoshConv(Double_t sign, Double_t tau, Double_t dgamma, Double_t sig, Double_t rtau, Double_t fsign) const
455 static Double_t root2(sqrt(2.)) ;
456 static Double_t root2pi(sqrt(2*atan2(0.,-1.))) ;
457 static Double_t rootpi(sqrt(atan2(0.,-1.)));
458 Double_t tau1 = 1/(1/tau-dgamma/2);
459 Double_t tau2 = 1/(1/tau+dgamma/2);
467 xp *= fsign ; // modified FMV,08/13/03
468 sign *= fsign ; // modified FMV,08/13/03
470 cFly=tau1*(exp(sig*sig/(2*tau1*tau1)-sign*xp/tau1)
471 *RooMath::erfc(sig/(root2*tau1)-sign*xp/(root2*sig))
472 +sign*exp(sig*sig/(2*rtau*rtau)+xp/rtau)
473 *RooMath::erfc(sig/(root2*rtau)+xp/(root2*sig)))/(2*(tau1+sign*rtau))
474 +tau2*(exp(sig*sig/(2*tau2*tau2)-sign*xp/tau2)
475 *RooMath::erfc(sig/(root2*tau2)-sign*xp/(root2*sig))
476 +sign*exp(sig*sig/(2*rtau*rtau)+xp/rtau)
477 *RooMath::erfc(sig/(root2*rtau)+xp/(root2*sig)))/(2*(tau2+sign*rtau));;
486 Double_t RooGExpModel::calcSinhConv(Double_t sign, Double_t sign1, Double_t sign2, Double_t tau, Double_t dgamma, Double_t sig, Double_t rtau, Double_t fsign) const
488 static Double_t root2(sqrt(2.)) ;
489 static Double_t root2pi(sqrt(2*atan2(0.,-1.))) ;
490 static Double_t rootpi(sqrt(atan2(0.,-1.)));
491 Double_t tau1 = 1/(1/tau-dgamma/2);
492 Double_t tau2 = 1/(1/tau+dgamma/2);
501 xp *= fsign ; // modified FMV,08/13/03
502 sign1 *= fsign ; // modified FMV,08/13/03
503 sign2 *= fsign ; // modified FMV,08/13/03
505 cFly=sign1*tau1*(exp(sig*sig/(2*tau1*tau1)-sign*xp/tau1)
506 *RooMath::erfc(sig/(root2*tau1)-sign*xp/(root2*sig))
507 +sign*exp(sig*sig/(2*rtau*rtau)+xp/rtau)
508 *RooMath::erfc(sig/(root2*rtau)+xp/(root2*sig)))/(2*(tau1+sign*rtau))
509 +sign2*tau2*(exp(sig*sig/(2*tau2*tau2)-sign*xp/tau2)
510 *RooMath::erfc(sig/(root2*tau2)-sign*xp/(root2*sig))
511 +sign*exp(sig*sig/(2*rtau*rtau)+xp/rtau)
512 *RooMath::erfc(sig/(root2*rtau)+xp/(root2*sig)))/(2*(tau2+sign*rtau));;
584 if (dGamma==0) basisType =
expBasis;
592 if (
verboseEval()>0) cout <<
"RooGExpModel::analyticalIntegral(" <<
GetName() <<
") 1st form" << endl ;
605 result = 0.5*
evalCerfInt(-fsign,rtau,fsign*umin,fsign*umax,c)/rtau ;
610 return result*ssfInt ;
618 if (
verboseEval()>0) cout <<
"RooGExpModel::analyticalIntegral(" <<
GetName() <<
") 2nd form" << endl ;
631 return result*ssfInt ;
637 if (
verboseEval()>0) cout <<
"RooGExpModel::analyticalIntegral(" <<
GetName() <<
") 4th form omega = "
638 << omega <<
", tau = " << tau << endl ;
641 if (wt==0)
return result ;
645 if (basisSign!=
Minus) result += -1*
calcSinConvNorm(+1,tau,omega,sig,rtau,fsign,rangeName).imag();
646 if (basisSign!=
Plus) result += -1*
calcSinConvNorm(-1,tau,omega,sig,rtau,fsign,rangeName).imag();
648 return result*ssfInt ;
654 <<
") 5th form omega = " << omega <<
", tau = " << tau << endl ;
661 if (basisSign!=
Plus) result +=
calcSinConvNorm(-1,tau,omega,sig,rtau,fsign,rangeName).real();
663 return result*ssfInt ;
671 <<
") 6th form dgamma = " << dgamma <<
", tau = " << tau << endl ;
690 <<
") 6th form dgamma = " << dgamma <<
", tau = " << tau << endl ;
735 std::complex<Double_t> eins(1,0);
736 std::complex<Double_t> k(1/tau,sign*omega);
737 std::complex<Double_t> term1 =
evalCerfInt(sign,-sign*omega*tau, tau, -sign*umin1, -sign*umax1, c1);
739 std::complex<Double_t> term2 = std::complex<Double_t>(
evalCerfInt(-fsign, rtau, fsign*umin2, fsign*umax2, c2)*fsign*
sign,0);
740 return (term1+term2)/(eins + k*fsign*sign*rtau) ;
769 if (
fabs(tau-rtau)<1e-10 &&
fabs(term1+term2)<1e-10) {
770 cout <<
"epsilon method" << endl ;
772 return calcSinConvNorm(sign,tau+epsilon,sig,rtau-epsilon,fsign,rangeName) ;
774 return (term1+term2)/(eins + k*fsign*sign*rtau) ;
784 std::complex<Double_t> diff;
786 diff = std::complex<Double_t>(2,0) ;
790 return std::complex<Double_t>(tau/(1.+wt*wt),0)*std::complex<Double_t>(1,wt)*diff;
802 if ((umin<-8 && umax>8)||(umax<-8 && umin>8)) {
822 std::complex<Double_t> z(swt*c,u+c);
823 std::complex<Double_t> zc(u+c,-swt*c);
824 std::complex<Double_t> zsq= z*z;
825 std::complex<Double_t>
v= -zsq - u*u;
836 if (
matchArgs(directVars,generateVars,
x))
return 1 ;
XYZVector ans(TestRotation const &t, XYZVector const &v_in)
Double_t logErfC(Double_t x) const
Approximation of the log of the complex error function.
std::complex< Double_t > evalCerfInt(Double_t sign, Double_t wt, Double_t tau, Double_t umin, Double_t umax, Double_t c) const
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...
void generateEvent(Int_t code)
Interface for generation of anan event using the algorithm corresponding to the specified code...
Double_t calcDecayConv(Double_t sign, Double_t tau, Double_t sig, Double_t rtau, Double_t fsign) const
Double_t evalCerfRe(Double_t u, Double_t c) const
virtual Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
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...
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
const RooAbsReal & arg() const
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
static std::complex< Double_t > evalCerfApprox(Double_t swt, Double_t u, Double_t c)
use the approximation: erf(z) = exp(-z*z)/(sqrt(pi)*z) to explicitly cancel the divergent exp(y*y) be...
std::complex< Double_t > calcSinConvNorm(Double_t sign, Double_t tau, Double_t omega, Double_t sig, Double_t rtau, Double_t fsign, const char *rangeName) const
old code (asymptotic normalization only) std::complex z(1/tau,sign*omega); return z*2/(omeg...
virtual Int_t basisCode(const char *name) const
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
virtual Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
Double_t getVal(const RooArgSet *set=0) const
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
RooRealVar & convVar() const
Return the convolution variable of the resolution model.
return fString CompareTo(((TObjString *) obj) ->fString)
virtual const char * GetName() const
Returns name of object.
static std::complex< double > erfc(const std::complex< double > z)
complex erfc function
static std::complex< Double_t > evalCerf(Double_t swt, Double_t u, Double_t c)
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
const RooFormulaVar & basis() const
double atan2(double, double)
virtual ~RooGExpModel()
Destructor.
virtual Double_t evaluate() const
std::complex< Double_t > calcSinConv(Double_t sign, Double_t sig, Double_t tau, Double_t omega, Double_t rtau, Double_t fsign) const
Double_t min(const char *rname=0) const
static std::complex< double > erf(const std::complex< double > z)
complex erf function
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Double_t max(const char *rname=0) const