JHUGen MELA  v2.4.1
Matrix element calculations as used in JHUGen. MELA is an important tool that was used for the Higgs boson discovery and for precise measurements of its structure and interactions. Please see the website https://spin.pha.jhu.edu/ and papers cited there for more details, and kindly cite those papers when using this code.
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
mela.m
m
Definition: mela.py:715
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