JHUGen MELA  JHUGen v7.5.6, MELA v2.4.2
Matrix element calculations as used in JHUGen.
RooSpinZero_3D_pp_VH.cc
Go to the documentation of this file.
1 #include "RooSpinZero_3D_pp_VH.h"
2 
3 
6  const char *name, const char *title,
7  modelMeasurables const& _measurables,
8  modelParameters const& _parameters,
9  modelCouplings const& _couplings,
10  Double_t _sqrts,
11  RooSpin::VdecayType _Vdecay1, RooSpin::VdecayType _Vdecay2,
12  TVar::VerbosityLevel verbosity_
13 ) : RooSpinZero(
14  name, title,
15  _measurables,
16  _parameters,
17  _couplings,
18  _Vdecay1, _Vdecay2,
19  verbosity_
20 ),
21 sqrts(_sqrts)
22 {}
23 
24 
26  const RooSpinZero_3D_pp_VH& other, const char* name
27 ) : RooSpinZero(other, name),
28 sqrts(other.sqrts)
29 {}
30 
31 
33  Double_t& A00term, Double_t& Appterm, Double_t& Ammterm,
34  Double_t& A00ppterm, Double_t& A00mmterm, Double_t& Appmmterm,
35  const Int_t code,
36  int VGammaVpmode1, int VGammaVpmode2
37 ) const{
38  const Double_t Pi = TMath::Pi();
39 
40  Double_t R1Val, R2Val;
41  calculateVffR1R2(R1Val, R2Val, VGammaVpmode1==1, VGammaVpmode2==1);
42  if (VGammaVpmode1==2 || VGammaVpmode2==2){
43  Double_t RVp1Val=0, RVp2Val=0;
44  calculateVprimeffR1R2(RVp1Val, RVp2Val);
45  if (VGammaVpmode1==2) R1Val=RVp1Val;
46  if (VGammaVpmode2==2) R1Val=RVp2Val;
47  }
48 
49  Double_t A00Re, A00Im, AppRe, AppIm, AmmRe, AmmIm;
50  calculateAmplitudes(A00Re, A00Im, AppRe, AppIm, AmmRe, AmmIm, VGammaVpmode1, VGammaVpmode2);
51 
52  Double_t A00 = A00Im*A00Im + A00Re*A00Re;
53  Double_t App = AppIm*AppIm + AppRe*AppRe;
54  Double_t Amm = AmmIm*AmmIm + AmmRe*AmmRe;
55 
56  Double_t phi00=atan2(A00Im, A00Re);
57  Double_t phipp=atan2(AppIm, AppRe)-phi00;
58  Double_t phimm=atan2(AmmIm, AmmRe)-phi00;
59 
60  Double_t A00_prefactor = 1.;
61  Double_t Amm_pp_prefactor = 1.;
62  Double_t A00mm_prefactor = 2.;
63  Double_t A00pp_prefactor = 2.;
64  Double_t Ammpp_prefactor = 2.;
65 
66  A00term = A00*A00_prefactor;
67  Appterm = App*Amm_pp_prefactor;
68  Ammterm = Amm*Amm_pp_prefactor;
69  A00ppterm = sqrt(A00*App)*A00pp_prefactor;
70  A00mmterm = sqrt(A00*Amm)*A00mm_prefactor;
71  Appmmterm = sqrt(Amm*App)*Ammpp_prefactor;
72 
73  if ((code % prime_h1)==0){
74  Double_t A00_h1int = 4./3.;
75  Double_t Ammpp_h1int = 8./3.;
76  Double_t A0m_h1int = Pi/2.*R1Val;
77  Double_t A0p_h1int = Pi/2.*R1Val;
78  Double_t Amp_h1int = 4./3.;
79 
80  if (A00term!=0) A00term *= A00_h1int;
81  if (Appterm!=0) Appterm *= Ammpp_h1int;
82  if (Ammterm!=0) Ammterm *= Ammpp_h1int;
83  if (A00ppterm!=0) A00ppterm *= A0p_h1int;
84  if (A00mmterm!=0) A00mmterm *= A0m_h1int;
85  if (Appmmterm!=0) Appmmterm *= Amp_h1int;
86  }
87  else{
88  if (A00term!=0) A00term *= (1. - pow(h1, 2));
89  if (Appterm!=0) Appterm *= (1. + pow(h1, 2) - 2.*h1*R1Val);
90  if (Ammterm!=0) Ammterm *= (1. + pow(h1, 2) + 2.*h1*R1Val);
91  if (A00ppterm!=0) A00ppterm *= sqrt(1. - pow(h1, 2))*(R1Val - h1);
92  if (A00mmterm!=0) A00mmterm *= sqrt(1. - pow(h1, 2))*(R1Val + h1);
93  if (Appmmterm!=0) Appmmterm *= (1. - pow(h1, 2));
94  }
95 
96  if ((code % prime_h2)==0){
97  Double_t A00_h2int = 4./3.;
98  Double_t Ammpp_h2int = 8./3.;
99  Double_t A0m_h2int = Pi/2.*R2Val;
100  Double_t A0p_h2int = Pi/2.*R2Val;
101  Double_t Amp_h2int = 4./3.;
102 
103  if (A00term!=0) A00term *= A00_h2int;
104  if (Appterm!=0) Appterm *= Ammpp_h2int;
105  if (Ammterm!=0) Ammterm *= Ammpp_h2int;
106  if (A00ppterm!=0) A00ppterm *= A0p_h2int;
107  if (A00mmterm!=0) A00mmterm *= A0m_h2int;
108  if (Appmmterm!=0) Appmmterm *= Amp_h2int;
109  }
110  else{
111  if (A00term!=0) A00term *= (1. - pow(h2, 2));
112  if (Appterm!=0) Appterm *= (1. + pow(h2, 2) + 2.*h2*R2Val);
113  if (Ammterm!=0) Ammterm *= (1. + pow(h2, 2) - 2.*h2*R2Val);
114  if (A00ppterm!=0) A00ppterm *= sqrt(1. - pow(h2, 2))*(R2Val + h2);
115  if (A00mmterm!=0) A00mmterm *= sqrt(1. - pow(h2, 2))*(R2Val - h2);
116  if (Appmmterm!=0) Appmmterm *= (1. - pow(h2, 2));
117  }
118 
119  if ((code % prime_hs)==0){
120  Double_t A_hsint = 2.;
121 
122  if (A00term!=0) A00term *= A_hsint;
123  if (Appterm!=0) Appterm *= A_hsint;
124  if (Ammterm!=0) Ammterm *= A_hsint;
125  if (A00ppterm!=0) A00ppterm *= A_hsint;
126  if (A00mmterm!=0) A00mmterm *= A_hsint;
127  if (Appmmterm!=0) Appmmterm *= A_hsint;
128  } // else *= 1
129 
130  if ((code % prime_Phi)==0){
131  Double_t A00mmpp_phiint = 2.*Pi;
132  Double_t A0p_phiint = 0;
133  Double_t A0m_phiint = 0;
134  Double_t Amp_phiint = 0;
135 
136  if (A00term!=0) A00term *= A00mmpp_phiint;
137  if (Appterm!=0) Appterm *= A00mmpp_phiint;
138  if (Ammterm!=0) Ammterm *= A00mmpp_phiint;
139  if (A00ppterm!=0) A00ppterm *= A0p_phiint;
140  if (A00mmterm!=0) A00mmterm *= A0m_phiint;
141  if (Appmmterm!=0) Appmmterm *= Amp_phiint;
142  }
143  else{
144  //if (A00term!=0) A00term *= 1.;
145  //if (Appterm!=0) Appterm *= 1.;
146  //if (Ammterm!=0) Ammterm *= 1.;
147  if (A00ppterm!=0) A00ppterm *= cos(Phi + phipp);
148  if (A00mmterm!=0) A00mmterm *= cos(Phi - phimm);
149  if (Appmmterm!=0) Appmmterm *= cos(2*Phi - phimm + phipp);
150  }
151 
152  if ((code % prime_Phi1)==0){
153  Double_t A_phi1int = 2.*Pi;
154 
155  if (A00term!=0) A00term *= A_phi1int;
156  if (Appterm!=0) Appterm *= A_phi1int;
157  if (Ammterm!=0) Ammterm *= A_phi1int;
158  if (A00ppterm!=0) A00ppterm *= A_phi1int;
159  if (A00mmterm!=0) A00mmterm *= A_phi1int;
160  if (Appmmterm!=0) Appmmterm *= A_phi1int;
161  }
162  // else *= 1
163 }
164 
165 
167  Double_t epsilon=1e-15;
168  Double_t m1_=m1; if (Vdecay1==RooSpin::kVdecayType_GammaOnshell) return epsilon;
169  Double_t m2_=m2; if (Vdecay2==RooSpin::kVdecayType_GammaOnshell) m2_=0;
170  if ((m12+m2_) > m1_ || (m2_ <= 0. && Vdecay2!=RooSpin::kVdecayType_GammaOnshell) || m1_ <= 0.) return epsilon;
171 
172  Int_t code = intCodeStart;
174 
175  Double_t betaValSq = (1.-(pow(m12-m2_, 2)/pow(m1_, 2)))*(1.-(pow(m12+m2_, 2)/pow(m1_, 2)));
176  if (betaValSq<=0.) return epsilon;
177  Double_t betaVal = sqrt(betaValSq);
178 
179  Double_t term1Coeff = 1;
180  Double_t term2Coeff = 1;
181  term1Coeff = pow(m1_*GeVunit, -2);
182  if (Vdecay2!=RooSpin::kVdecayType_GammaOnshell) term2Coeff = 2.*m2_*GeVunit;
183  Double_t plumi = partonicLuminosity(m1_, Y, sqrts);
184 
185  Double_t value = 0;
186  for (int VGammaVpmode1=0; VGammaVpmode1<=2; VGammaVpmode1++){
187  for (int VGammaVpmode2=0; VGammaVpmode2<=2; VGammaVpmode2++){
188  if (!(
190  &&
191  (VGammaVpmode2==1 || Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
192  )
193  ||
194  VGammaVpmode1==1
195  ||
196  (VGammaVpmode1==2 && VGammaVpmode2==1)
197  ||
198  !computeNeededAmplitude(VGammaVpmode1, VGammaVpmode2)
199  ) continue;
200  Double_t val_A00=0, val_App=0, val_Amm=0, val_A0p=0, val_A0m=0, val_Amp=0;
201  evaluatePolarizationTerms(val_A00, val_App, val_Amm, val_A0p, val_A0m, val_Amp, code, VGammaVpmode1, VGammaVpmode2);
202  value += val_A00 + val_App + val_Amm + val_A0p + val_A0m + val_Amp;
203  }
204  }
205  value *= betaVal*term1Coeff*term2Coeff*plumi;
206  return value;
207 }
208 
209 Int_t RooSpinZero_3D_pp_VH::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const{
210  Int_t code = intCodeStart;
211  if (checkFundamentalType(h1)){ if (matchArgs(allVars, analVars, h1)) code *= prime_h1; }
212  if (checkFundamentalType(h2)){ if (matchArgs(allVars, analVars, h2) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_h2; }
213  if (checkFundamentalType(hs)){ if (matchArgs(allVars, analVars, hs)) code *= prime_hs; }
214  if (checkFundamentalType(Phi)){ if (matchArgs(allVars, analVars, Phi) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_Phi; }
215  if (checkFundamentalType(Phi1)){ if (matchArgs(allVars, analVars, Phi1) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_Phi1; }
216  if (code==1) code=0;
217  return code;
218 }
219 Double_t RooSpinZero_3D_pp_VH::analyticalIntegral(Int_t code, const char* /*rangeName*/) const{
220  Double_t epsilon=1e-10;
221  Double_t m1_=m1; if (Vdecay1==RooSpin::kVdecayType_GammaOnshell) return epsilon;
222  Double_t m2_=m2; if (Vdecay2==RooSpin::kVdecayType_GammaOnshell) m2_=0;
223  if ((m12+m2_) > m1_ || (m2_ <= 0. && Vdecay2!=RooSpin::kVdecayType_GammaOnshell) || m1_ <= 0.) return epsilon;
224 
225  Double_t betaValSq = (1.-(pow(m12-m2_, 2)/pow(m1_, 2)))*(1.-(pow(m12+m2_, 2)/pow(m1_, 2)));
226  if (betaValSq<=0.) return epsilon;
227  Double_t betaVal = sqrt(betaValSq);
228 
229  Double_t term1Coeff = 1;
230  Double_t term2Coeff = 1;
231  term1Coeff = pow(m1_*GeVunit, -2);
232  if (Vdecay2!=RooSpin::kVdecayType_GammaOnshell) term2Coeff = 2.*m2_*GeVunit;
233  Double_t plumi = partonicLuminosity(m1_, Y, sqrts);
234 
235  Double_t value = 0;
236  for (int VGammaVpmode1=0; VGammaVpmode1<=2; VGammaVpmode1++){
237  for (int VGammaVpmode2=0; VGammaVpmode2<=2; VGammaVpmode2++){
238  if (!(
239  (VGammaVpmode1==0 && VGammaVpmode2==0 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
240  ||
241  (VGammaVpmode1==0 && VGammaVpmode2==1 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell)
242  ||
243  (VGammaVpmode1==1 && VGammaVpmode2==0 && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
244  ||
245  (VGammaVpmode1==1 && VGammaVpmode2==1)
246  ||
247  (VGammaVpmode1==0 && VGammaVpmode2==2 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
248  ||
249  (VGammaVpmode1==2 && VGammaVpmode2==0 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
250  ||
251  (VGammaVpmode1==2 && VGammaVpmode2==2 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
252  )) continue;
253  Double_t val_A00=0, val_App=0, val_Amm=0, val_A0p=0, val_A0m=0, val_Amp=0;
254  evaluatePolarizationTerms(val_A00, val_App, val_Amm, val_A0p, val_A0m, val_Amp, code, VGammaVpmode1, VGammaVpmode2);
255  value += val_A00 + val_App + val_Amm + val_A0p + val_A0m + val_Amp;
256  }
257  }
258  value *= betaVal*term1Coeff*term2Coeff*plumi;
259  return value;
260 }
261 
262 Double_t RooSpinZero_3D_pp_VH::partonicLuminosity(Double_t mVal, Double_t YVal, Double_t sqrtsVal) const{
263  Double_t Q = mVal;
264  Double_t xa0 = mVal/sqrtsVal; // -> E E <- => xa0=xb0=p/Ebeam (Ebeam = sqrts/2)
265  Double_t xa = exp(YVal)*xa0; // -> Ea Eb <- => xa=Ea/Ebeam
266  Double_t xb = exp(-YVal)*xa0; // -> Ea Eb <- => xb=Eb/Ebeam
267  Double_t prefactor = 1.;
268  if (
269  (
270  (mVal <= 600. && fabs(YVal) > 20.*pow(mVal, -0.32))
271  ||
272  (mVal > 600. && fabs(YVal) > 21.*pow(mVal, -0.34))
274  ){
275  xa = xa0;
276  xb = xa0;
277  prefactor = 1e-5;
278  }
279 
280  Double_t weightu = 0.5;
281  Double_t weightd = 0.5;
282  Double_t weightc = 1.0;
283  Double_t weights = 1.0;
284  Double_t weightb = 1.0;
285 
286 
287  // PDF parameters
288  // up params
289  Double_t u0par0 = 0.03134; Double_t u0par1 =-2.068e-05; Double_t u0par2 = 1.283e-08;
290  Double_t u1par0 = 0.9; Double_t u1par1 =-0.0004307; Double_t u1par2 = 2.458e-07;
291  Double_t u2par0 =-0.1369; Double_t u2par1 = 0.003423; Double_t u2par2 =-2.155e-06;
292  Double_t u3par0 =-0.4013; Double_t u3par1 =-0.0002574; Double_t u3par2 = 1.561e-07;
293  Double_t u4par0 = 0.5782; Double_t u4par1 =-0.004728; Double_t u4par2 = 2.906e-06;
294  Double_t ubar0par0 = 0.02856; Double_t ubar0par1 =-2.112e-05; Double_t ubar0par2 = 1.272e-08;
295  Double_t ubar1par0 =-0.06822; Double_t ubar1par1 = 3.172e-05; Double_t ubar1par2 =-2.008e-08;
296  Double_t ubar2par0 = 0.1967; Double_t ubar2par1 =-0.000118; Double_t ubar2par2 = 6.871e-08;
297  Double_t ubar3par0 =-0.2251; Double_t ubar3par1 = 0.0001295; Double_t ubar3par2 =-7.181e-08;
298  Double_t ubar4par0 =-0.4068; Double_t ubar4par1 =-0.0002956; Double_t ubar4par2 = 1.783e-07;
299  Double_t ubar5par0 =-2.251; Double_t ubar5par1 =-0.0001699; Double_t ubar5par2 = 1.492e-07;
300  // down params
301  Double_t d0par0 = 0.03278; Double_t d0par1 =-2.915e-05; Double_t d0par2 = 1.809e-08;
302  Double_t d1par0 = 0.479; Double_t d1par1 =-0.0002559; Double_t d1par2 = 1.557e-07;
303  Double_t d2par0 =-0.5972; Double_t d2par1 = 0.0003118; Double_t d2par2 =-1.905e-07;
304  Double_t d3par0 =-0.3892; Double_t d3par1 =-0.000317; Double_t d3par2 = 1.944e-07;
305  Double_t d4par0 = 0.5007; Double_t d4par1 =-0.001665; Double_t d4par2 = 9.895e-07;
306  Double_t dbar0par0 = 0.02328; Double_t dbar0par1 =-1.367e-05; Double_t dbar0par2 = 8.246e-09;
307  Double_t dbar1par0 = 0.09422; Double_t dbar1par1 =-0.0001019; Double_t dbar1par2 = 6.375e-08;
308  Double_t dbar2par0 =-0.5296; Double_t dbar2par1 = 0.000466; Double_t dbar2par2 =-2.896e-07;
309  Double_t dbar3par0 = 0.5354; Double_t dbar3par1 =-0.0004404; Double_t dbar3par2 = 2.728e-07;
310  Double_t dbar4par0 =-0.4386; Double_t dbar4par1 =-0.0002605; Double_t dbar4par2 = 1.582e-07;
311  Double_t dbar5par0 =-1.289; Double_t dbar5par1 =-0.001618; Double_t dbar5par2 = 9.601e-07;
312  // charm, strange, bottom params
313  Double_t c0par0 = 0.01829; Double_t c0par1 =-6.93e-06; Double_t c0par2 = 3.796e-09;
314  Double_t c1par0 = 0.03081; Double_t c1par1 = 4.325e-05; Double_t c1par2 =-3.95e-08;
315  Double_t c2par0 = 0.5398; Double_t c2par1 =-4.284e-05; Double_t c2par2 =-1.362e-08;
316  Double_t c3par0 =-0.5986; Double_t c3par1 = 0.002565; Double_t c3par2 =-1.937e-06;
317  Double_t c4par0 =-0.4534; Double_t c4par1 =-0.0002329; Double_t c4par2 = 1.343e-07;
318  Double_t c5par0 =-8.657; Double_t c5par1 =-0.005157; Double_t c5par2 = 3.68e-06;
319  Double_t s0par0 = 0.01312; Double_t s0par1 =-3.743e-06; Double_t s0par2 = 2.076e-09;
320  Double_t s1par0 =-0.001416; Double_t s1par1 =-7.649e-06; Double_t s1par2 = 4.757e-09;
321  Double_t s2par0 = 0.2864; Double_t s2par1 =-6.693e-05; Double_t s2par2 = 3.566e-08;
322  Double_t s3par0 =-0.4857; Double_t s3par1 =-0.000253; Double_t s3par2 = 1.541e-07;
323  Double_t s4par0 =-10.33; Double_t s4par1 =-0.001601; Double_t s4par2 = 9.718e-07;
324  Double_t b0par0 = 0.005934; Double_t b0par1 = 2.516e-06; Double_t b0par2 =-1.828e-09;
325  Double_t b1par0 =-0.003063; Double_t b1par1 =-6.761e-06; Double_t b1par2 = 4.298e-09;
326  Double_t b2par0 = 0.1174; Double_t b2par1 = 3.752e-05; Double_t b2par2 =-2.863e-08;
327  Double_t b3par0 =-0.5549; Double_t b3par1 =-0.0002205; Double_t b3par2 = 1.334e-07;
328  Double_t b4par0 =-10.18; Double_t b4par1 =-0.001136; Double_t b4par2 = 6.931e-07;
329 
330  // PDF definition
331  Double_t up0 = u0par0 + u0par1*Q + u0par2*Q*Q;
332  Double_t up1 = u1par0 + u1par1*Q + u1par2*Q*Q;
333  Double_t up2 = u2par0 + u2par1*Q + u2par2*Q*Q;
334  Double_t up3 = u3par0 + u3par1*Q + u3par2*Q*Q;
335  Double_t up4 = u4par0 + u4par1*Q + u4par2*Q*Q;
336  Double_t antiup0 = ubar0par0 + ubar0par1*Q + ubar0par2*Q*Q;
337  Double_t antiup1 = ubar1par0 + ubar1par1*Q + ubar1par2*Q*Q;
338  Double_t antiup2 = ubar2par0 + ubar2par1*Q + ubar2par2*Q*Q;
339  Double_t antiup3 = ubar3par0 + ubar3par1*Q + ubar3par2*Q*Q;
340  Double_t antiup4 = ubar4par0 + ubar4par1*Q + ubar4par2*Q*Q;
341  Double_t antiup5 = ubar5par0 + ubar5par1*Q + ubar5par2*Q*Q;
342  Double_t down0 = d0par0 + d0par1*Q + d0par2*Q*Q;
343  Double_t down1 = d1par0 + d1par1*Q + d1par2*Q*Q;
344  Double_t down2 = d2par0 + d2par1*Q + d2par2*Q*Q;
345  Double_t down3 = d3par0 + d3par1*Q + d3par2*Q*Q;
346  Double_t down4 = d4par0 + d4par1*Q + d4par2*Q*Q;
347  Double_t antidown0 = dbar0par0 + dbar0par1*Q + dbar0par2*Q*Q;
348  Double_t antidown1 = dbar1par0 + dbar1par1*Q + dbar1par2*Q*Q;
349  Double_t antidown2 = dbar2par0 + dbar2par1*Q + dbar2par2*Q*Q;
350  Double_t antidown3 = dbar3par0 + dbar3par1*Q + dbar3par2*Q*Q;
351  Double_t antidown4 = dbar4par0 + dbar4par1*Q + dbar4par2*Q*Q;
352  Double_t antidown5 = dbar5par0 + dbar5par1*Q + dbar5par2*Q*Q;
353  Double_t charm0 = c0par0 + c0par1*Q + c0par2*Q*Q;
354  Double_t charm1 = c1par0 + c1par1*Q + c1par2*Q*Q;
355  Double_t charm2 = c2par0 + c2par1*Q + c2par2*Q*Q;
356  Double_t charm3 = c3par0 + c3par1*Q + c3par2*Q*Q;
357  Double_t charm4 = c4par0 + c4par1*Q + c4par2*Q*Q;
358  Double_t charm5 = c5par0 + c5par1*Q + c5par2*Q*Q;
359  Double_t strange0 = s0par0 + s0par1*Q + s0par2*Q*Q;
360  Double_t strange1 = s1par0 + s1par1*Q + s1par2*Q*Q;
361  Double_t strange2 = s2par0 + s2par1*Q + s2par2*Q*Q;
362  Double_t strange3 = s3par0 + s3par1*Q + s3par2*Q*Q;
363  Double_t strange4 = s4par0 + s4par1*Q + s4par2*Q*Q;
364  Double_t bottom0 = b0par0 + b0par1*Q + b0par2*Q*Q;
365  Double_t bottom1 = b1par0 + b1par1*Q + b1par2*Q*Q;
366  Double_t bottom2 = b2par0 + b2par1*Q + b2par2*Q*Q;
367  Double_t bottom3 = b3par0 + b3par1*Q + b3par2*Q*Q;
368  Double_t bottom4 = b4par0 + b4par1*Q + b4par2*Q*Q;
369 
370  Double_t FuncAu1 = (up0+up1*xa+up2*pow(xa, 2))*pow((1-xa), 4)*pow(xa, up3)*exp(1.0+up4*xa);
371  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);
372  Double_t FuncAu2 = (up0+up1*xb+up2*pow(xb, 2))*pow((1-xb), 4)*pow(xb, up3)*exp(1.0+up4*xb);
373  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);
374 
375  Double_t FuncAd1 = (down0+down1*xa+down2*pow(xa, 2))*pow((1-xa), 4)*pow(xa, down3)*exp(1.0+down4*xa);
376  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);
377  Double_t FuncAd2 = (down0+down1*xb+down2*pow(xb, 2))*pow((1-xb), 4)*pow(xb, down3)*exp(1.0+down4*xb);
378  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);
379 
380  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);
381  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);
382  Double_t Funcsa = (strange0+strange1*xa+strange2*pow(xa, 2))*pow((1-xa), 4)*pow(xa, strange3)*exp(1.0+strange4*xa);
383  Double_t Funcsb = (strange0+strange1*xb+strange2*pow(xb, 2))*pow((1-xb), 4)*pow(xb, strange3)*exp(1.0+strange4*xb);
384  Double_t Funcba = (bottom0+bottom1*xa+bottom2*pow(xa, 2))*pow((1-xa), 4)*pow(xa, bottom3)*exp(1.0+bottom4*xa);
385  Double_t Funcbb = (bottom0+bottom1*xb+bottom2*pow(xb, 2))*pow((1-xb), 4)*pow(xb, bottom3)*exp(1.0+bottom4*xb);
386 
387  Double_t totSec = 0;
388 
389  if (Vdecay1!=RooSpin::kVdecayType_Wany){ // ZH or gammaH
390  Double_t FuncABu = FuncAu1/xa*FuncBu1/xb+FuncAu2/xa*FuncBu2/xb;
391  Double_t FuncABd = FuncAd1/xa*FuncBd1/xb+FuncAd2/xa*FuncBd2/xb;
392  Double_t FuncABc = Funcsa*Funcsb/xa/xb;
393  Double_t FuncABs = Funcca*Funccb/xa/xb;
394  Double_t FuncABb = Funcba*Funcbb/xa/xb;
395 
396  totSec = 2.*prefactor*mVal*(
397  (FuncABu)*weightu
398  +(FuncABd)*weightd
399  +(FuncABc)*weightc
400  +(FuncABs)*weights
401  +(FuncABb)*weightb
402  );
403  }
404  else{ // WH
405  Double_t FuncAB_udbar = FuncAu1/xa*FuncBd1/xb + FuncAu2/xa*FuncBd2/xb;
406  Double_t FuncAB_csbar = Funcca*Funcsb/xa/xb;
407  Double_t FuncAB_dubar = FuncAd1/xa*FuncBu1/xb + FuncAd2/xa*FuncBu2/xb;
408  Double_t FuncAB_scbar = Funccb*Funcsa/xa/xb;
409 
410  totSec = 2.*prefactor*mVal*(
411  (FuncAB_udbar)*weightu
412  +(FuncAB_dubar)*weightd
413  +(FuncAB_csbar)*weightc
414  +(FuncAB_scbar)*weights
415  );
416  }
417 
418  if (totSec<=0) totSec = 1e-5;
419  return totSec;
420 }
421 
422 
423 
value
pymela::gHIGGS_KAPPA value("gHIGGS_KAPPA_TILDE", pymela::gHIGGS_KAPPA_TILDE) .value("SIZE_HQQ"
RooSpin::Phi
RooRealProxy Phi
Definition: RooSpin.h:110
TVar::VerbosityLevel
VerbosityLevel
Definition: TVar.hh:47
RooSpin::prime_hs
@ prime_hs
Definition: RooSpin.h:41
RooSpin::m12
RooRealProxy m12
Definition: RooSpin.h:113
RooSpinZero_3D_pp_VH.h
RooSpin::m1
RooRealProxy m1
Definition: RooSpin.h:111
RooSpinZero::computeNeededAmplitude
virtual Bool_t computeNeededAmplitude(int VGammaVpmode1, int VGammaVpmode2) const final
Definition: RooSpinZero.cc:583
RooSpin::calculateVprimeffR1R2
virtual void calculateVprimeffR1R2(Double_t &R1Val, Double_t &R2Val) const
Definition: RooSpin.cc:343
RooSpin::prime_Phi
@ prime_Phi
Definition: RooSpin.h:42
RooSpin::checkFundamentalType
virtual Bool_t checkFundamentalType(const RooRealProxy &proxy) const
Definition: RooSpin.cc:392
RooSpin::prime_h2
@ prime_h2
Definition: RooSpin.h:40
RooSpin::Y
RooRealProxy Y
Definition: RooSpin.h:116
RooSpin::Vdecay2
RooSpin::VdecayType Vdecay2
Definition: RooSpin.h:136
RooSpinZero
Definition: RooSpinZero.h:7
RooSpinZero_3D_pp_VH::analyticalIntegral
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Definition: RooSpinZero_3D_pp_VH.cc:219
RooSpinZero_3D_pp_VH
Definition: RooSpinZero_3D_pp_VH.h:9
RooSpin::intCodeStart
Int_t intCodeStart
Definition: RooSpin.h:138
RooSpin::kVdecayType_GammaOnshell
@ kVdecayType_GammaOnshell
Definition: RooSpin.h:30
RooSpinZero_3D_pp_VH::sqrts
Double_t sqrts
Definition: RooSpinZero_3D_pp_VH.h:12
hto_aux_hbb::xb
real *8 xb
Definition: CALLING_cpHTO.f:226
RooSpinZero_3D_pp_VH::partonicLuminosity
Double_t partonicLuminosity(Double_t mVal, Double_t YVal, Double_t sqrts) const
Definition: RooSpinZero_3D_pp_VH.cc:262
RooSpinZero::modelCouplings
Definition: RooSpinZero.h:10
RooSpinZero_3D_pp_VH::evaluate
Double_t evaluate() const
Definition: RooSpinZero_3D_pp_VH.cc:166
RooSpin::m2
RooRealProxy m2
Definition: RooSpin.h:112
RooSpin::h2
RooRealProxy h2
Definition: RooSpin.h:109
RooSpinZero::calculateAmplitudes
virtual void calculateAmplitudes(Double_t &A00Re, Double_t &A00Im, Double_t &AppRe, Double_t &AppIm, Double_t &AmmRe, Double_t &AmmIm, int VGammaVpmode1=0, int VGammaVpmode2=0) const
Definition: RooSpinZero.cc:477
RooSpinZero_3D_pp_VH::evaluatePolarizationTerms
void evaluatePolarizationTerms(Double_t &A00term, Double_t &Appterm, Double_t &Ammterm, Double_t &A00ppterm, Double_t &A00mmterm, Double_t &Appmmterm, const Int_t code, int VGammaVpmode1=0, int VGammaVpmode2=0) const
Definition: RooSpinZero_3D_pp_VH.cc:32
RooSpin::hs
RooRealProxy hs
Definition: RooSpin.h:114
RooSpinZero_3D_pp_VH::getAnalyticalIntegral
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Definition: RooSpinZero_3D_pp_VH.cc:209
RooSpin::prime_h1
@ prime_h1
Definition: RooSpin.h:39
RooSpin::GeVunit
static constexpr Double_t GeVunit
Definition: RooSpin.h:141
Q
double Q[11]
Definition: TMCFM.hh:121
RooSpin::prime_Phi1
@ prime_Phi1
Definition: RooSpin.h:43
RooSpin::Vdecay1
RooSpin::VdecayType Vdecay1
Definition: RooSpin.h:135
RooSpin::kVdecayType_Wany
@ kVdecayType_Wany
Definition: RooSpin.h:29
RooSpin::modelMeasurables
Definition: RooSpin.h:50
sqrts
double sqrts
Definition: TMCFM.hh:290
RooSpin::h1
RooRealProxy h1
Definition: RooSpin.h:108
RooSpin::modelParameters
Definition: RooSpin.h:61
RooSpinZero_3D_pp_VH::RooSpinZero_3D_pp_VH
RooSpinZero_3D_pp_VH()
Definition: RooSpinZero_3D_pp_VH.cc:4
RooSpin::calculateVffR1R2
virtual void calculateVffR1R2(Double_t &R1Val, Double_t &R2Val, bool isGammaV1=false, bool isGammaV2=false) const
Definition: RooSpin.cc:297
RooSpin::VdecayType
VdecayType
Definition: RooSpin.h:28
RooSpin::Phi1
RooRealProxy Phi1
Definition: RooSpin.h:115