JHUGen MELA  JHUGen v7.5.6, MELA v2.4.2
Matrix element calculations as used in JHUGen.
RooRapidityBkg.cc
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * *
4  * This code was autogenerated by RooClassFactory *
5  *****************************************************************************/
6 
7 // Your description goes here...
8 
9 #include <math.h>
10 #include "RooRapidityBkg.h"
11 #include "Riostream.h"
12 #include "RooAbsReal.h"
13 #include "RooAbsCategory.h"
14 //#include "RooComplex.h"
15 #include "RooMath.h"
16 #include "TMath.h"
17 
18 //ClassImp(RooRapidityBkg)
19 
20 RooRapidityBkg::RooRapidityBkg(const char *name, const char *title,
21  RooAbsReal& _Y,
22  RooAbsReal& _m,
23  RooAbsReal& _sqrtS) :
24 RooAbsPdf(name,title),
25 Y("Y","Y",this,_Y),
26 m("m","m",this,_m),
27 sqrtS("sqrtS","sqrtS",this,_sqrtS)
28 {
29 }
30 
31 
32 RooRapidityBkg::RooRapidityBkg(const RooRapidityBkg& other, const char* name) :
33 RooAbsPdf(other,name),
34 Y("Y",this,other.Y),
35 m("m",this,other.m),
36 sqrtS("sqrtS",this,other.sqrtS)
37 {
38 }
39 
40 
41 
42 Double_t RooRapidityBkg::evaluate() const
43 {
44  Double_t s0 = sqrtS*sqrtS;
45  Double_t s = m*m;
46  Double_t Q = m;
47  Double_t xa = exp(Y)*sqrt(s/s0);
48  Double_t xb = exp(-Y)*sqrt(s/s0);
49 
50 
51  Double_t weightu = 0.5;
52  Double_t weightd = 0.5;
53  Double_t weightc = 1.0;
54  Double_t weights = 1.0;
55  Double_t weightb = 1.0;
56 
57 
58 
59  // PDF parameters
60  // up params
61  Double_t u0par0 = 0.03134; Double_t u0par1 =-2.068e-05; Double_t u0par2 = 1.283e-08;
62  Double_t u1par0 = 0.9; Double_t u1par1 =-0.0004307; Double_t u1par2 = 2.458e-07;
63  Double_t u2par0 =-0.1369; Double_t u2par1 = 0.003423; Double_t u2par2 =-2.155e-06;
64  Double_t u3par0 =-0.4013; Double_t u3par1 =-0.0002574; Double_t u3par2 = 1.561e-07;
65  Double_t u4par0 = 0.5782; Double_t u4par1 =-0.004728; Double_t u4par2 = 2.906e-06;
66  Double_t ubar0par0 = 0.02856; Double_t ubar0par1 =-2.112e-05; Double_t ubar0par2 = 1.272e-08;
67  Double_t ubar1par0 =-0.06822; Double_t ubar1par1 = 3.172e-05; Double_t ubar1par2 =-2.008e-08;
68  Double_t ubar2par0 = 0.1967; Double_t ubar2par1 =-0.000118; Double_t ubar2par2 = 6.871e-08;
69  Double_t ubar3par0 =-0.2251; Double_t ubar3par1 = 0.0001295; Double_t ubar3par2 =-7.181e-08;
70  Double_t ubar4par0 =-0.4068; Double_t ubar4par1 =-0.0002956; Double_t ubar4par2 = 1.783e-07;
71  Double_t ubar5par0 =-2.251; Double_t ubar5par1 =-0.0001699; Double_t ubar5par2 = 1.492e-07;
72  // down params
73  Double_t d0par0 = 0.03278; Double_t d0par1 =-2.915e-05; Double_t d0par2 = 1.809e-08;
74  Double_t d1par0 = 0.479; Double_t d1par1 =-0.0002559; Double_t d1par2 = 1.557e-07;
75  Double_t d2par0 =-0.5972; Double_t d2par1 = 0.0003118; Double_t d2par2 =-1.905e-07;
76  Double_t d3par0 =-0.3892; Double_t d3par1 =-0.000317; Double_t d3par2 = 1.944e-07;
77  Double_t d4par0 = 0.5007; Double_t d4par1 =-0.001665; Double_t d4par2 = 9.895e-07;
78  Double_t dbar0par0 = 0.02328; Double_t dbar0par1 =-1.367e-05; Double_t dbar0par2 = 8.246e-09;
79  Double_t dbar1par0 = 0.09422; Double_t dbar1par1 =-0.0001019; Double_t dbar1par2 = 6.375e-08;
80  Double_t dbar2par0 =-0.5296; Double_t dbar2par1 = 0.000466; Double_t dbar2par2 =-2.896e-07;
81  Double_t dbar3par0 = 0.5354; Double_t dbar3par1 =-0.0004404; Double_t dbar3par2 = 2.728e-07;
82  Double_t dbar4par0 =-0.4386; Double_t dbar4par1 =-0.0002605; Double_t dbar4par2 = 1.582e-07;
83  Double_t dbar5par0 =-1.289; Double_t dbar5par1 =-0.001618; Double_t dbar5par2 = 9.601e-07;
84  // charm, strange, bottom params
85  Double_t c0par0 = 0.01829; Double_t c0par1 =-6.93e-06; Double_t c0par2 = 3.796e-09;
86  Double_t c1par0 = 0.03081; Double_t c1par1 = 4.325e-05; Double_t c1par2 =-3.95e-08;
87  Double_t c2par0 = 0.5398; Double_t c2par1 =-4.284e-05; Double_t c2par2 =-1.362e-08;
88  Double_t c3par0 =-0.5986; Double_t c3par1 = 0.002565; Double_t c3par2 =-1.937e-06;
89  Double_t c4par0 =-0.4534; Double_t c4par1 =-0.0002329; Double_t c4par2 = 1.343e-07;
90  Double_t c5par0 =-8.657; Double_t c5par1 =-0.005157; Double_t c5par2 = 3.68e-06;
91  Double_t s0par0 = 0.01312; Double_t s0par1 =-3.743e-06; Double_t s0par2 = 2.076e-09;
92  Double_t s1par0 =-0.001416; Double_t s1par1 =-7.649e-06; Double_t s1par2 = 4.757e-09;
93  Double_t s2par0 = 0.2864; Double_t s2par1 =-6.693e-05; Double_t s2par2 = 3.566e-08;
94  Double_t s3par0 =-0.4857; Double_t s3par1 =-0.000253; Double_t s3par2 = 1.541e-07;
95  Double_t s4par0 =-10.33; Double_t s4par1 =-0.001601; Double_t s4par2 = 9.718e-07;
96  Double_t b0par0 = 0.005934; Double_t b0par1 = 2.516e-06; Double_t b0par2 =-1.828e-09;
97  Double_t b1par0 =-0.003063; Double_t b1par1 =-6.761e-06; Double_t b1par2 = 4.298e-09;
98  Double_t b2par0 = 0.1174; Double_t b2par1 = 3.752e-05; Double_t b2par2 =-2.863e-08;
99  Double_t b3par0 =-0.5549; Double_t b3par1 =-0.0002205; Double_t b3par2 = 1.334e-07;
100  Double_t b4par0 =-10.18; Double_t b4par1 =-0.001136; Double_t b4par2 = 6.931e-07;
101 
102 
103  // PDF definition
104  Double_t up0 = u0par0 + u0par1*Q + u0par2*Q*Q;
105  Double_t up1 = u1par0 + u1par1*Q + u1par2*Q*Q;
106  Double_t up2 = u2par0 + u2par1*Q + u2par2*Q*Q;
107  Double_t up3 = u3par0 + u3par1*Q + u3par2*Q*Q;
108  Double_t up4 = u4par0 + u4par1*Q + u4par2*Q*Q;
109  Double_t antiup0 = ubar0par0 + ubar0par1*Q + ubar0par2*Q*Q;
110  Double_t antiup1 = ubar1par0 + ubar1par1*Q + ubar1par2*Q*Q;
111  Double_t antiup2 = ubar2par0 + ubar2par1*Q + ubar2par2*Q*Q;
112  Double_t antiup3 = ubar3par0 + ubar3par1*Q + ubar3par2*Q*Q;
113  Double_t antiup4 = ubar4par0 + ubar4par1*Q + ubar4par2*Q*Q;
114  Double_t antiup5 = ubar5par0 + ubar5par1*Q + ubar5par2*Q*Q;
115  Double_t down0 = d0par0 + d0par1*Q + d0par2*Q*Q;
116  Double_t down1 = d1par0 + d1par1*Q + d1par2*Q*Q;
117  Double_t down2 = d2par0 + d2par1*Q + d2par2*Q*Q;
118  Double_t down3 = d3par0 + d3par1*Q + d3par2*Q*Q;
119  Double_t down4 = d4par0 + d4par1*Q + d4par2*Q*Q;
120  Double_t antidown0 = dbar0par0 + dbar0par1*Q + dbar0par2*Q*Q;
121  Double_t antidown1 = dbar1par0 + dbar1par1*Q + dbar1par2*Q*Q;
122  Double_t antidown2 = dbar2par0 + dbar2par1*Q + dbar2par2*Q*Q;
123  Double_t antidown3 = dbar3par0 + dbar3par1*Q + dbar3par2*Q*Q;
124  Double_t antidown4 = dbar4par0 + dbar4par1*Q + dbar4par2*Q*Q;
125  Double_t antidown5 = dbar5par0 + dbar5par1*Q + dbar5par2*Q*Q;
126  Double_t charm0 = c0par0 + c0par1*Q + c0par2*Q*Q;
127  Double_t charm1 = c1par0 + c1par1*Q + c1par2*Q*Q;
128  Double_t charm2 = c2par0 + c2par1*Q + c2par2*Q*Q;
129  Double_t charm3 = c3par0 + c3par1*Q + c3par2*Q*Q;
130  Double_t charm4 = c4par0 + c4par1*Q + c4par2*Q*Q;
131  Double_t charm5 = c5par0 + c5par1*Q + c5par2*Q*Q;
132  Double_t strange0 = s0par0 + s0par1*Q + s0par2*Q*Q;
133  Double_t strange1 = s1par0 + s1par1*Q + s1par2*Q*Q;
134  Double_t strange2 = s2par0 + s2par1*Q + s2par2*Q*Q;
135  Double_t strange3 = s3par0 + s3par1*Q + s3par2*Q*Q;
136  Double_t strange4 = s4par0 + s4par1*Q + s4par2*Q*Q;
137  Double_t bottom0 = b0par0 + b0par1*Q + b0par2*Q*Q;
138  Double_t bottom1 = b1par0 + b1par1*Q + b1par2*Q*Q;
139  Double_t bottom2 = b2par0 + b2par1*Q + b2par2*Q*Q;
140  Double_t bottom3 = b3par0 + b3par1*Q + b3par2*Q*Q;
141  Double_t bottom4 = b4par0 + b4par1*Q + b4par2*Q*Q;
142 
143  Double_t FuncAu1 = (up0+up1*xa+up2*pow(xa,2))*pow((1-xa),4)*pow(xa,up3)*exp(1.0+up4*xa);
144  Double_t FuncBu1 = (antiup0+antiup1*xb+antiup2*pow(xb,2)+antiup3*pow(xb,3))*pow((1-xb),4)*pow(xb,antiup4)*exp(1.0+antiup5*xb);
145  Double_t FuncAu2 = (up0+up1*xb+up2*pow(xb,2))*pow((1-xb),4)*pow(xb,up3)*exp(1.0+up4*xb);
146  Double_t FuncBu2 = (antiup0+antiup1*xa+antiup2*pow(xa,2)+antiup3*pow(xa,3))*pow((1-xa),4)*pow(xa,antiup4)*exp(1.0+antiup5*xa);
147  Double_t FuncABu = FuncAu1/xa*FuncBu1/xb+FuncAu2/xa*FuncBu2/xb;
148 
149 
150  Double_t FuncAd1 = (down0+down1*xa+down2*pow(xa,2))*pow((1-xa),4)*pow(xa,down3)*exp(1.0+down4*xa);
151  Double_t FuncBd1 = (antidown0+antidown1*xb+antidown2*pow(xb,2)+antidown3*pow(xb,3))*pow((1-xb),4)*pow(xb,antidown4)*exp(1.0+antidown5*xb);
152  Double_t FuncAd2 = (down0+down1*xb+down2*pow(xb,2))*pow((1-xb),4)*pow(xb,down3)*exp(1.0+down4*xb);
153  Double_t FuncBd2 = (antidown0+antidown1*xa+antidown2*pow(xa,2)+antidown3*pow(xa,3))*pow((1-xa),4)*pow(xa,antidown4)*exp(1.0+antidown5*xa);
154  Double_t FuncABd = FuncAd1/xa*FuncBd1/xb+FuncAd2/xa*FuncBd2/xb;
155 
156  Double_t Funcca = (charm0+charm1*xa+charm2*pow(xa,2)+charm3*pow(xa,3))*pow((1-xa),4)*pow(xa,charm4)*exp(1.0+charm5*xa);
157  Double_t Funccb = (charm0+charm1*xb+charm2*pow(xb,2)+charm3*pow(xb,3))*pow((1-xb),4)*pow(xb,charm4)*exp(1.0+charm5*xb);
158  Double_t Funcsa = (strange0+strange1*xa+strange2*pow(xa,2))*pow((1-xa),4)*pow(xa,strange3)*exp(1.0+strange4*xa);
159  Double_t Funcsb = (strange0+strange1*xb+strange2*pow(xb,2))*pow((1-xb),4)*pow(xb,strange3)*exp(1.0+strange4*xb);
160  Double_t Funcba = (bottom0+bottom1*xa+bottom2*pow(xa,2))*pow((1-xa),4)*pow(xa,bottom3)*exp(1.0+bottom4*xa);
161  Double_t Funcbb = (bottom0+bottom1*xb+bottom2*pow(xb,2))*pow((1-xb),4)*pow(xb,bottom3)*exp(1.0+bottom4*xb);
162  Double_t FuncABc = Funcsa*Funcsb/xa/xb;
163  Double_t FuncABs = Funcca*Funccb/xa/xb;
164  Double_t FuncABb = Funcba*Funcbb/xa/xb;
165 
166 
167  Double_t totSec = 2*m*(
168  (FuncABu)*weightu
169  +(FuncABd)*weightd
170  +(FuncABc)*weightc
171  +(FuncABs)*weights
172  +(FuncABb)*weightb
173  );
174 
175  if(( m <= 600. && TMath::Abs(Y) > 20*pow(m,-0.32)) || ( m > 600. && TMath::Abs(Y) > 21*pow(m,-0.34)))
176  {
177  //Find totSec when mZZ, Y=0
178  Double_t xa0 = sqrt(s/s0); //at Y=0 xa=xb
179 
180  //up
181  //if xa=xb then FuncAu1=FuncAu2 and FuncBu1=FuncBu2
182  FuncAu1 = (up0+up1*xa0+up2*pow(xa0,2))*pow((1-xa0),4)*pow(xa0,up3)*exp(1.0+up4*xa0);
183  FuncBu1 = (antiup0+antiup1*xa0+antiup2*pow(xa0,2)+antiup3*pow(xa0,3))*pow((1-xa0),4)*pow(xa0,antiup4)*exp(1.0+antiup5*xa0);
184  FuncABu = 2*(FuncAu1/xa0*FuncBu1/xa0);
185 
186  //down
187  //if xa=xb then FuncAd1=FuncAd2 and FuncBd1=FuncBd2
188  FuncAd1 = (down0+down1*xa0+down2*pow(xa0,2))*pow((1-xa0),4)*pow(xa0,down3)*exp(1.0+down4*xa0);
189  FuncBd1 = (antidown0+antidown1*xa0+antidown2*pow(xa0,2)+antidown3*pow(xa0,3))*pow((1-xa0),4)*pow(xa0,antidown4)*exp(1.0+antidown5*xa0);
190  FuncABd = 2*(FuncAd1/xa0*FuncBd1/xa0);
191 
192  //sea
193  Funcca = (charm0+charm1*xa0+charm2*pow(xa0,2)+charm3*pow(xa0,3))*pow((1-xa0),4)*pow(xa0,charm4)*exp(1.0+charm5*xa0); //Funcca=Funccb
194  Funcsa = (strange0+strange1*xa0+strange2*pow(xa0,2))*pow((1-xa0),4)*pow(xa0,strange3)*exp(1.0+strange4*xa0); //Funcsa=Funcsb
195  Funcba = (bottom0+bottom1*xa0+bottom2*pow(xa0,2))*pow((1-xa0),4)*pow(xa0,bottom3)*exp(1.0+bottom4*xa0); //Funcba=Funcbb
196  FuncABc = Funcsa*Funcsa/xa0/xa0;
197  FuncABs = Funcca*Funcca/xa0/xa0;
198  FuncABb = Funcba*Funcba/xa0/xa0;
199  Double_t totSec0 = 2*m*(
200  (FuncABu)*weightu
201  +(FuncABd)*weightd
202  +(FuncABc)*weightc
203  +(FuncABs)*weights
204  +(FuncABb)*weightb
205  );
206  totSec = 1.e-5*totSec0;
207  }
208 
209  if (totSec <= 0.) totSec = 0.00001;
210  return totSec;
211 
212 }
213 
214 
215 
216 
217 
218 
RooRapidityBkg::Y
RooRealProxy Y
Definition: RooRapidityBkg.h:32
RooRapidityBkg::sqrtS
RooRealProxy sqrtS
Definition: RooRapidityBkg.h:34
RooRapidityBkg::evaluate
Double_t evaluate() const
Definition: RooRapidityBkg.cc:42
RooRapidityBkg::RooRapidityBkg
RooRapidityBkg()
Definition: RooRapidityBkg.h:22
hto_aux_hbb::xb
real *8 xb
Definition: CALLING_cpHTO.f:226
Q
double Q[11]
Definition: TMCFM.hh:121
RooRapidityBkg::m
RooRealProxy m
Definition: RooRapidityBkg.h:33
RooRapidityBkg
Definition: RooRapidityBkg.h:16
RooRapidityBkg.h