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.
RooSpinZero_5D_VH.cc
Go to the documentation of this file.
1 #include "RooSpinZero_5D_VH.h"
2 
3 
6  const char *name, const char *title,
7  modelMeasurables const& _measurables,
8  modelParameters const& _parameters,
9  modelCouplings const& _couplings,
10  RooSpin::VdecayType _Vdecay1, RooSpin::VdecayType _Vdecay2,
11  TVar::VerbosityLevel verbosity_
12 ) : RooSpinZero(
13  name, title,
14  _measurables,
15  _parameters,
16  _couplings,
17  _Vdecay1, _Vdecay2,
18  verbosity_
19 )
20 {}
21 
22 
24  const RooSpinZero_5D_VH& other, const char* name
25 ) : RooSpinZero(other, name)
26 {}
27 
29  Double_t& A00term, Double_t& Appterm, Double_t& Ammterm,
30  Double_t& A00ppterm, Double_t& A00mmterm, Double_t& Appmmterm,
31  const Int_t code,
32  int VGammaVpmode1, int VGammaVpmode2
33 ) const{
34  const Double_t Pi = TMath::Pi();
35 
36  Double_t R1Val, R2Val;
37  calculateVffR1R2(R1Val, R2Val, VGammaVpmode1==1, VGammaVpmode2==1);
38  if (VGammaVpmode1==2 || VGammaVpmode2==2){
39  Double_t RVp1Val=0, RVp2Val=0;
40  calculateVprimeffR1R2(RVp1Val, RVp2Val);
41  if (VGammaVpmode1==2) R1Val=RVp1Val;
42  if (VGammaVpmode2==2) R1Val=RVp2Val;
43  }
44 
45  Double_t A00Re, A00Im, AppRe, AppIm, AmmRe, AmmIm;
46  calculateAmplitudes(A00Re, A00Im, AppRe, AppIm, AmmRe, AmmIm, VGammaVpmode1, VGammaVpmode2);
47 
48  Double_t A00 = A00Im*A00Im + A00Re*A00Re;
49  Double_t App = AppIm*AppIm + AppRe*AppRe;
50  Double_t Amm = AmmIm*AmmIm + AmmRe*AmmRe;
51 
52  Double_t phi00=atan2(A00Im, A00Re);
53  Double_t phipp=atan2(AppIm, AppRe)-phi00;
54  Double_t phimm=atan2(AmmIm, AmmRe)-phi00;
55 
56  Double_t A00_prefactor = 1.;
57  Double_t Amm_pp_prefactor = 1.;
58  Double_t A00mm_prefactor = 2.;
59  Double_t A00pp_prefactor = 2.;
60  Double_t Ammpp_prefactor = 2.;
61 
62  A00term = A00*A00_prefactor;
63  Appterm = App*Amm_pp_prefactor;
64  Ammterm = Amm*Amm_pp_prefactor;
65  A00ppterm = sqrt(A00*App)*A00pp_prefactor;
66  A00mmterm = sqrt(A00*Amm)*A00mm_prefactor;
67  Appmmterm = sqrt(Amm*App)*Ammpp_prefactor;
68 
69  if ((code % prime_h1)==0){
70  Double_t A00_h1int = 4./3.;
71  Double_t Ammpp_h1int = 8./3.;
72  Double_t A0m_h1int = Pi/2.*R1Val;
73  Double_t A0p_h1int = Pi/2.*R1Val;
74  Double_t Amp_h1int = 4./3.;
75 
76  if (A00term!=0) A00term *= A00_h1int;
77  if (Appterm!=0) Appterm *= Ammpp_h1int;
78  if (Ammterm!=0) Ammterm *= Ammpp_h1int;
79  if (A00ppterm!=0) A00ppterm *= A0p_h1int;
80  if (A00mmterm!=0) A00mmterm *= A0m_h1int;
81  if (Appmmterm!=0) Appmmterm *= Amp_h1int;
82  }
83  else{
84  if (A00term!=0) A00term *= (1. - pow(h1, 2));
85  if (Appterm!=0) Appterm *= (1. + pow(h1, 2) - 2.*h1*R1Val);
86  if (Ammterm!=0) Ammterm *= (1. + pow(h1, 2) + 2.*h1*R1Val);
87  if (A00ppterm!=0) A00ppterm *= sqrt(1. - pow(h1, 2))*(R1Val - h1);
88  if (A00mmterm!=0) A00mmterm *= sqrt(1. - pow(h1, 2))*(R1Val + h1);
89  if (Appmmterm!=0) Appmmterm *= (1. - pow(h1, 2));
90  }
91 
92  if ((code % prime_h2)==0){
93  Double_t A00_h2int = 4./3.;
94  Double_t Ammpp_h2int = 8./3.;
95  Double_t A0m_h2int = Pi/2.*R2Val;
96  Double_t A0p_h2int = Pi/2.*R2Val;
97  Double_t Amp_h2int = 4./3.;
98 
99  if (A00term!=0) A00term *= A00_h2int;
100  if (Appterm!=0) Appterm *= Ammpp_h2int;
101  if (Ammterm!=0) Ammterm *= Ammpp_h2int;
102  if (A00ppterm!=0) A00ppterm *= A0p_h2int;
103  if (A00mmterm!=0) A00mmterm *= A0m_h2int;
104  if (Appmmterm!=0) Appmmterm *= Amp_h2int;
105  }
106  else{
107  if (A00term!=0) A00term *= (1. - pow(h2, 2));
108  if (Appterm!=0) Appterm *= (1. + pow(h2, 2) + 2.*h2*R2Val);
109  if (Ammterm!=0) Ammterm *= (1. + pow(h2, 2) - 2.*h2*R2Val);
110  if (A00ppterm!=0) A00ppterm *= sqrt(1. - pow(h2, 2))*(R2Val + h2);
111  if (A00mmterm!=0) A00mmterm *= sqrt(1. - pow(h2, 2))*(R2Val - h2);
112  if (Appmmterm!=0) Appmmterm *= (1. - pow(h2, 2));
113  }
114 
115  if ((code % prime_hs)==0){
116  Double_t A_hsint = 2.;
117 
118  if (A00term!=0) A00term *= A_hsint;
119  if (Appterm!=0) Appterm *= A_hsint;
120  if (Ammterm!=0) Ammterm *= A_hsint;
121  if (A00ppterm!=0) A00ppterm *= A_hsint;
122  if (A00mmterm!=0) A00mmterm *= A_hsint;
123  if (Appmmterm!=0) Appmmterm *= A_hsint;
124  } // else *= 1
125 
126  if ((code % prime_Phi)==0){
127  Double_t A00mmpp_phiint = 2.*Pi;
128  Double_t A0p_phiint = 0;
129  Double_t A0m_phiint = 0;
130  Double_t Amp_phiint = 0;
131 
132  if (A00term!=0) A00term *= A00mmpp_phiint;
133  if (Appterm!=0) Appterm *= A00mmpp_phiint;
134  if (Ammterm!=0) Ammterm *= A00mmpp_phiint;
135  if (A00ppterm!=0) A00ppterm *= A0p_phiint;
136  if (A00mmterm!=0) A00mmterm *= A0m_phiint;
137  if (Appmmterm!=0) Appmmterm *= Amp_phiint;
138  }
139  else{
140  //if (A00term!=0) A00term *= 1.;
141  //if (Appterm!=0) Appterm *= 1.;
142  //if (Ammterm!=0) Ammterm *= 1.;
143  if (A00ppterm!=0) A00ppterm *= cos(Phi + phipp);
144  if (A00mmterm!=0) A00mmterm *= cos(Phi - phimm);
145  if (Appmmterm!=0) Appmmterm *= cos(2*Phi - phimm + phipp);
146  }
147 
148  if ((code % prime_Phi1)==0){
149  Double_t A_phi1int = 2.*Pi;
150 
151  if (A00term!=0) A00term *= A_phi1int;
152  if (Appterm!=0) Appterm *= A_phi1int;
153  if (Ammterm!=0) Ammterm *= A_phi1int;
154  if (A00ppterm!=0) A00ppterm *= A_phi1int;
155  if (A00mmterm!=0) A00mmterm *= A_phi1int;
156  if (Appmmterm!=0) Appmmterm *= A_phi1int;
157  }
158  // else *= 1
159 }
160 
162  Double_t epsilon=1e-15;
163  Double_t m1_=m1; if (Vdecay1==RooSpin::kVdecayType_GammaOnshell) return epsilon;
164  Double_t m2_=m2; if (Vdecay2==RooSpin::kVdecayType_GammaOnshell) m2_=0;
165  if ((m12+m2_) > m1_ || (m2_ <= 0. && Vdecay2!=RooSpin::kVdecayType_GammaOnshell) || m1_ <= 0.) return epsilon;
166 
167  Int_t code = intCodeStart;
169 
170  Double_t betaValSq = (1.-(pow(m12-m2_, 2)/pow(m1_, 2)))*(1.-(pow(m12+m2_, 2)/pow(m1_, 2)));
171  if (betaValSq<=0.) return epsilon;
172  Double_t betaVal = sqrt(betaValSq);
173 
174  Double_t term1Coeff = 1;
175  Double_t term2Coeff = 1;
176  term1Coeff = pow(m1_*GeVunit, -2);
177  if (Vdecay2!=RooSpin::kVdecayType_GammaOnshell) term2Coeff = 2.*m2_*GeVunit;
178 
179  Double_t value = 0;
180  for (int VGammaVpmode1=0; VGammaVpmode1<=2; VGammaVpmode1++){
181  for (int VGammaVpmode2=0; VGammaVpmode2<=2; VGammaVpmode2++){
182  if (!(
184  &&
185  (VGammaVpmode2==1 || Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
186  )
187  ||
188  VGammaVpmode1==1
189  ||
190  (VGammaVpmode1==2 && VGammaVpmode2==1)
191  ||
192  !computeNeededAmplitude(VGammaVpmode1, VGammaVpmode2)
193  ) continue;
194  Double_t val_A00=0, val_App=0, val_Amm=0, val_A0p=0, val_A0m=0, val_Amp=0;
195  evaluatePolarizationTerms(val_A00, val_App, val_Amm, val_A0p, val_A0m, val_Amp, code, VGammaVpmode1, VGammaVpmode2);
196  value += val_A00 + val_App + val_Amm + val_A0p + val_A0m + val_Amp;
197  }
198  }
199  value *= betaVal*term1Coeff*term2Coeff;
200  return value;
201 }
202 
203 Int_t RooSpinZero_5D_VH::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const{
204  Int_t code = intCodeStart;
205  if (checkFundamentalType(h1)){ if (matchArgs(allVars, analVars, h1)) code *= prime_h1; }
206  if (checkFundamentalType(h2)){ if (matchArgs(allVars, analVars, h2) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_h2; }
207  if (checkFundamentalType(hs)){ if (matchArgs(allVars, analVars, hs)) code *= prime_hs; }
208  if (checkFundamentalType(Phi)){ if (matchArgs(allVars, analVars, Phi) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_Phi; }
209  if (checkFundamentalType(Phi1)){ if (matchArgs(allVars, analVars, Phi1) || Vdecay2==RooSpin::kVdecayType_GammaOnshell) code *= prime_Phi1; }
210  if (code==1) code=0;
211  return code;
212 }
213 Double_t RooSpinZero_5D_VH::analyticalIntegral(Int_t code, const char* /*rangeName*/) const{
214  Double_t epsilon=1e-10;
215  Double_t m1_=m1; if (Vdecay1==RooSpin::kVdecayType_GammaOnshell) return epsilon;
216  Double_t m2_=m2; if (Vdecay2==RooSpin::kVdecayType_GammaOnshell) m2_=0;
217  if ((m12+m2_) > m1_ || (m2_ <= 0. && Vdecay2!=RooSpin::kVdecayType_GammaOnshell) || m1_ <= 0.) return epsilon;
218 
219  Double_t betaValSq = (1.-(pow(m12-m2_, 2)/pow(m1_, 2)))*(1.-(pow(m12+m2_, 2)/pow(m1_, 2)));
220  if (betaValSq<=0.) return epsilon;
221  Double_t betaVal = sqrt(betaValSq);
222 
223  Double_t term1Coeff = 1;
224  Double_t term2Coeff = 1;
225  term1Coeff = pow(m1_*GeVunit, -2);
226  if (Vdecay2!=RooSpin::kVdecayType_GammaOnshell) term2Coeff = 2.*m2_*GeVunit;
227 
228  Double_t value = 0;
229  for (int VGammaVpmode1=0; VGammaVpmode1<=2; VGammaVpmode1++){
230  for (int VGammaVpmode2=0; VGammaVpmode2<=2; VGammaVpmode2++){
231  if (!(
232  (VGammaVpmode1==0 && VGammaVpmode2==0 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
233  ||
234  (VGammaVpmode1==0 && VGammaVpmode2==1 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell)
235  ||
236  (VGammaVpmode1==1 && VGammaVpmode2==0 && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
237  ||
238  (VGammaVpmode1==1 && VGammaVpmode2==1)
239  ||
240  (VGammaVpmode1==0 && VGammaVpmode2==2 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
241  ||
242  (VGammaVpmode1==2 && VGammaVpmode2==0 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
243  ||
244  (VGammaVpmode1==2 && VGammaVpmode2==2 && Vdecay1!=RooSpin::kVdecayType_GammaOnshell && Vdecay2!=RooSpin::kVdecayType_GammaOnshell)
245  )) continue;
246  Double_t val_A00=0, val_App=0, val_Amm=0, val_A0p=0, val_A0m=0, val_Amp=0;
247  evaluatePolarizationTerms(val_A00, val_App, val_Amm, val_A0p, val_A0m, val_Amp, code, VGammaVpmode1, VGammaVpmode2);
248  value += val_A00 + val_App + val_Amm + val_A0p + val_A0m + val_Amp;
249  }
250  }
251  value *= betaVal*term1Coeff*term2Coeff;
252  return value;
253 }
254 
255 
256 
257 
RooSpinZero_5D_VH.h
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
RooSpinZero_5D_VH::RooSpinZero_5D_VH
RooSpinZero_5D_VH()
Definition: RooSpinZero_5D_VH.cc:4
RooSpin::m12
RooRealProxy m12
Definition: RooSpin.h:113
RooSpin::prime_h1
@ prime_h1
Definition: RooSpin.h:39
RooSpin::prime_Phi
@ prime_Phi
Definition: RooSpin.h:42
RooSpin::prime_Phi1
@ prime_Phi1
Definition: RooSpin.h:43
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::checkFundamentalType
virtual Bool_t checkFundamentalType(const RooRealProxy &proxy) const
Definition: RooSpin.cc:392
RooSpin::Vdecay2
RooSpin::VdecayType Vdecay2
Definition: RooSpin.h:136
RooSpinZero
Definition: RooSpinZero.h:7
RooSpinZero_5D_VH::evaluate
Double_t evaluate() const
Definition: RooSpinZero_5D_VH.cc:161
RooSpin::intCodeStart
Int_t intCodeStart
Definition: RooSpin.h:138
RooSpin::kVdecayType_GammaOnshell
@ kVdecayType_GammaOnshell
Definition: RooSpin.h:30
RooSpinZero_5D_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_5D_VH.cc:28
RooSpin::prime_h2
@ prime_h2
Definition: RooSpin.h:40
RooSpinZero::modelCouplings
Definition: RooSpinZero.h:10
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_5D_VH::analyticalIntegral
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Definition: RooSpinZero_5D_VH.cc:213
RooSpinZero_5D_VH
Definition: RooSpinZero_5D_VH.h:7
RooSpin::hs
RooRealProxy hs
Definition: RooSpin.h:114
RooSpin::GeVunit
static constexpr Double_t GeVunit
Definition: RooSpin.h:141
RooSpin::Vdecay1
RooSpin::VdecayType Vdecay1
Definition: RooSpin.h:135
RooSpin::modelMeasurables
Definition: RooSpin.h:50
RooSpinZero_5D_VH::getAnalyticalIntegral
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Definition: RooSpinZero_5D_VH.cc:203
RooSpin::h1
RooRealProxy h1
Definition: RooSpin.h:108
RooSpin::modelParameters
Definition: RooSpin.h:61
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::prime_hs
@ prime_hs
Definition: RooSpin.h:41
RooSpin::Phi1
RooRealProxy Phi1
Definition: RooSpin.h:115