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.
|
Go to the documentation of this file.
66 #include <sys/types.h>
78 #include "RooMsgService.h"
89 using namespace RooFit;
99 melaRandomNumber(35797),
101 myVerbosity_(verbosity_),
112 melaRandomNumber(35797),
113 LHCsqrts(other.LHCsqrts),
114 myVerbosity_(other.myVerbosity_),
163 remove(
"ffwarn.dat");
167 remove(
"process.DAT");
168 remove(
"Pdfdata/cteq6l1.tbl");
169 remove(
"Pdfdata/cteq6l.tbl");
170 remove(
"Pdfdata/NNPDF30_lo_as_0130.LHgrid");
180 const double maxSqrts = 8.;
188 const string mcfmWarning = MELAPKGPATH +
"data/ffwarn.dat"; symlink(mcfmWarning.c_str(),
"ffwarn.dat");
189 const string mcfm_brsm_o = MELAPKGPATH +
"data/br.sm1"; symlink(mcfm_brsm_o.c_str(),
"br.sm1");
190 const string mcfm_brsm_t = MELAPKGPATH +
"data/br.sm2"; symlink(mcfm_brsm_t.c_str(),
"br.sm2");
191 const string mcfmInput1 = MELAPKGPATH +
"data/input.DAT"; symlink(mcfmInput1.c_str(),
"input.DAT");
192 const string mcfmInput2 = MELAPKGPATH +
"data/process.DAT"; symlink(mcfmInput2.c_str(),
"process.DAT");
194 mkdir(
"Pdfdata", S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);
195 const string mcfmInput3 = MELAPKGPATH +
"data/Pdfdata/cteq6l1.tbl"; symlink(mcfmInput3.c_str(),
"Pdfdata/cteq6l1.tbl");
196 const string mcfmInput4 = MELAPKGPATH +
"data/Pdfdata/cteq6l.tbl"; symlink(mcfmInput4.c_str(),
"Pdfdata/cteq6l.tbl");
200 mzz_rrv =
new RooRealVar(
"mzz",
"m_{ZZ}", mh_, 0., 1000.);
201 z1mass_rrv =
new RooRealVar(
"z1mass",
"m_{Z1}", 0., 160.);
202 z2mass_rrv =
new RooRealVar(
"z2mass",
"m_{Z2}", 0., 200.);
203 costhetastar_rrv =
new RooRealVar(
"costhetastar",
"cos#theta^{*}", -1., 1.);
204 costheta1_rrv =
new RooRealVar(
"costheta1",
"cos#theta_{1}", -1., 1.);
205 costheta2_rrv =
new RooRealVar(
"costheta2",
"cos#theta_{2}", -1., 1.);
206 phi_rrv=
new RooRealVar(
"phi",
"#Phi", -TMath::Pi(), TMath::Pi());
207 phi1_rrv=
new RooRealVar(
"phi1",
"#Phi_{1}", -TMath::Pi(), TMath::Pi());
208 Y_rrv =
new RooRealVar(
"Yzz",
"#Y_{ZZ}", 0, -4, 4);
209 upFrac_rrv =
new RooRealVar(
"upFrac",
"fraction up-quarks", .5, 0., 1.);
226 qqZZmodel =
new RooqqZZ_JHU_ZgammaZZ_fast(
"qqZZmodel",
"qqZZmodel", *
z1mass_rrv, *
z2mass_rrv, *
costheta1_rrv, *
costheta2_rrv, *
phi_rrv, *
costhetastar_rrv, *
phi1_rrv, *
mzz_rrv, *
upFrac_rrv);
229 const string path_HiggsWidthFile = MELAPKGPATH +
"data/HiggsTotalWidth_YR3.txt";
231 const string path_nnpdf = MELAPKGPATH +
"data/Pdfdata/NNPDF30_lo_as_0130.LHgrid";
232 char path_nnpdf_c[] =
"Pdfdata/NNPDF30_lo_as_0130.LHgrid";
235 symlink(path_nnpdf.c_str(), path_nnpdf_c);
252 RooMsgService::instance().getStream(1).removeTopic(NumIntegration);
253 RooMsgService::instance().setStreamStatus(1, kFALSE);
254 RooMsgService::instance().setStreamStatus(0, kFALSE);
258 if (superMELA_LHCsqrts > maxSqrts) superMELA_LHCsqrts = maxSqrts;
261 sprintf(cardpath,
"data/CombinationInputs/SM_inputs_%dTeV/inputs_4mu.txt", superMELA_LHCsqrts);
262 string cardfile = MELAPKGPATH + cardpath;
270 if (superDijetSqrts<13.) superDijetSqrts=13.;
279 vector<string> logolines;
280 logolines.push_back(
"MELA (Matrix Element Likelihood Approach)");
281 logolines.push_back(
"");
282 logolines.push_back(
"Data analysis and Monte Carlo weights package");
283 logolines.push_back(
"for analyses of resonances produced at pp, ppbar, and e+e- colliders, featuring:");
284 logolines.push_back(
"");
285 logolines.push_back(
"* JHUGenMELA *");
286 logolines.push_back(
"Signal calculations based on analytical pdf.s, and JHU Generator (JHUGen) matrix elements");
287 logolines.push_back(
"(See JHUGen credits below)");
288 logolines.push_back(
"");
289 logolines.push_back(
"* MCFM *");
290 logolines.push_back(
"Signal, background, and interference calculations, modified based on JHUGen matrix elements");
291 logolines.push_back(
"(See MCFM credits below)");
292 logolines.push_back(
"");
293 logolines.push_back(
"For more details: http://spin.pha.jhu.edu");
294 logolines.push_back(
"");
295 size_t maxlinesize = 0;
296 for (
auto const&
l:logolines) maxlinesize = std::max(maxlinesize,
l.length());
298 unsigned int iline=0;
299 for (
auto const&
l:logolines){
433 for (
int im=0; im<2; im++){
458 for (
int im=0; im<2; im++){
468 for (
int im=0; im<2; im++){
474 for (
int im=0; im<2; im++){
485 for (
int im=0; im<2; im++){
497 for (
int im=0; im<2; im++){
502 for (
int im=0; im<2; im++){
524 std::vector<TLorentzVector>
Mela::calculate4Momentum(
double Mx,
double M1,
double M2,
double theta,
double theta1,
double theta2,
double Phi1,
double Phi){
571 qH=0; m1=0; m2=0; costheta1=0; costheta2=0; Phi=0; costhetastar=0; Phi1=0;
575 TLorentzVector
const nullVector(0, 0, 0, 0);
578 simple_event_record mela_event;
579 mela_event.AssociationCode=partIncCode;
583 if (
daughters.size()<2 ||
daughters.size()>4 || mela_event.intermediateVid.size()!=2){
584 if (
myVerbosity_>=
TVar::ERROR)
MELAerr <<
"Mela::computeDecayAngles: Number of daughters " <<
daughters.size() <<
" or number of intermediate Vs " << mela_event.intermediateVid <<
" not supported!" << endl;
599 costhetastar, costheta1, costheta2, Phi, Phi1,
607 if (!std::isfinite(costhetastar)) costhetastar=0;
608 if (!std::isfinite(costheta1)) costheta1=0;
609 if (!std::isfinite(costheta2)) costheta2=0;
610 if (!std::isfinite(Phi)) Phi=0;
611 if (!std::isfinite(Phi1)) Phi1=0;
632 Q2V1=0; Q2V2=0; costheta1=0; costheta2=0; Phi=0; costhetastar=0; Phi1=0;
636 TLorentzVector
const nullVector(0, 0, 0, 0);
638 int nRequested_AssociatedJets=2;
640 simple_event_record mela_event;
641 mela_event.AssociationCode=partIncCode;
642 mela_event.nRequested_AssociatedJets=nRequested_AssociatedJets;
648 if ((
int)aparts.size()!=nRequested_AssociatedJets){
if (
myVerbosity_>=
TVar::ERROR)
MELAerr <<
"Mela::computeVBFAngles: Number of associated particles is not 2!" << endl;
return; }
654 for (
auto it=
daughters.cbegin()+1; it!=
daughters.cend(); it++){ firstPart.second = firstPart.second + it->second; }
664 MELAout <<
"Mela::computeVBFAngles: Giving the following particles to TUtil::computeVBFAngles:" << endl;
665 for (
unsigned int i=0; i<std::min(
daughters.size(), (SimpleParticleCollection_t::size_type) 4); i++)
MELAout <<
daughters.at(i) << endl;
666 for (
unsigned int i=0; i<std::min(aparts.size(), (SimpleParticleCollection_t::size_type) 2); i++)
MELAout << aparts.at(i) << endl;
667 for (
unsigned int i=0; i<std::min(
mothers.size(), (SimpleParticleCollection_t::size_type) 2); i++)
MELAout <<
mothers.at(i) << endl;
671 costhetastar, costheta1, costheta2, Phi, Phi1, Q2V1, Q2V2,
676 aparts.at(0).second, aparts.at(0).first,
677 aparts.at(1).second, aparts.at(1).first,
683 if (!std::isfinite(costhetastar)) costhetastar=0;
684 if (!std::isfinite(costheta1)) costheta1=0;
685 if (!std::isfinite(costheta2)) costheta2=0;
686 if (!std::isfinite(Phi)) Phi=0;
687 if (!std::isfinite(Phi1)) Phi1=0;
688 if (!std::isfinite(Q2V1)) Q2V1=0;
689 if (!std::isfinite(Q2V2)) Q2V2=0;
692 <<
"Mela::computeVBFAngles: (Q2_1, Q2_2, h1, h2, Phi, hs, Phi1) = "
693 << Q2V1 <<
", " << Q2V2 <<
", "
694 << costheta1 <<
", " << costheta2 <<
", " << Phi <<
", "
695 << costhetastar <<
", " << Phi1 << endl;
705 float& costheta1_real,
float& costheta1_imag,
706 float& costheta2_real,
float& costheta2_imag,
714 Q2V1=0; Q2V2=0; costheta1_real=0; costheta2_real=0; costheta1_imag=0; costheta2_imag=0; Phi=0; costhetastar=0; Phi1=0;
718 TLorentzVector
const nullVector(0, 0, 0, 0);
720 int nRequested_AssociatedJets=2;
722 simple_event_record mela_event;
723 mela_event.AssociationCode=partIncCode;
724 mela_event.nRequested_AssociatedJets=nRequested_AssociatedJets;
730 if ((
int) aparts.size()!=nRequested_AssociatedJets){
if (
myVerbosity_>=
TVar::ERROR)
MELAerr <<
"Mela::computeVBFAngles_ComplexBoost: Number of associated particles is not 2!" << endl;
return; }
736 for (
auto it=
daughters.cbegin()+1; it!=
daughters.cend(); it++){ firstPart.second = firstPart.second + it->second; }
746 costhetastar, costheta1_real, costheta1_imag, costheta2_real, costheta2_imag, Phi, Phi1, Q2V1, Q2V2,
751 aparts.at(0).second, aparts.at(0).first,
752 aparts.at(1).second, aparts.at(1).first,
758 if (!std::isfinite(costhetastar)) costhetastar=0;
759 if (!std::isfinite(costheta1_real)) costheta1_real=0;
760 if (!std::isfinite(costheta2_real)) costheta2_real=0;
761 if (!std::isfinite(costheta1_imag)) costheta1_imag=0;
762 if (!std::isfinite(costheta2_imag)) costheta2_imag=0;
763 if (!std::isfinite(Phi)) Phi=0;
764 if (!std::isfinite(Phi1)) Phi1=0;
765 if (!std::isfinite(Q2V1)) Q2V1=0;
766 if (!std::isfinite(Q2V2)) Q2V2=0;
769 <<
", " << costheta1_real <<
" + " << costheta1_imag <<
"i, "
770 << costheta2_real <<
" + " << costheta2_imag <<
"i, " << Phi <<
", "
771 << costhetastar <<
", " << Phi1 << endl;
773 <<
"Mela::computeVBFAngles_ComplexBoost: (Q2_1, Q2_2, h1, h2, Phi, hs, Phi1) = "
774 << Q2V1 <<
", " << Q2V2 <<
", "
775 << costheta1_real <<
" + " << costheta1_imag <<
"i, "
776 << costheta2_real <<
" + " << costheta2_imag <<
"i, "
777 << Phi <<
", " << costhetastar <<
", " << Phi1 << endl;
779 else if (
myVerbosity_>=
TVar::DEBUG)
MELAerr <<
"Mela::computeVBFAngles_ComplexBoost: No possible melaCand in TEvtProb to compute angles." << endl;
798 mVstar = 0; mV = 0; costheta1=0; costheta2=0; Phi=0; costhetastar=0; Phi1=0;
802 TLorentzVector
const nullVector(0, 0, 0, 0);
809 int nRequested_AssociatedJets=0;
810 int nRequested_AssociatedLeptons=0;
811 int nRequested_AssociatedPhotons=0;
812 int AssociationVCompatibility=0;
816 nRequested_AssociatedJets=2;
820 nRequested_AssociatedLeptons=2;
824 nRequested_AssociatedPhotons=1;
829 simple_event_record mela_event;
830 mela_event.AssociationCode=partIncCode;
831 mela_event.AssociationVCompatibility=AssociationVCompatibility;
832 mela_event.nRequested_AssociatedJets=nRequested_AssociatedJets;
833 mela_event.nRequested_AssociatedLeptons=nRequested_AssociatedLeptons;
834 mela_event.nRequested_AssociatedPhotons=nRequested_AssociatedPhotons;
842 MELAerr <<
"Mela::computeVHAngles: Number of associated particles (" << aparts.size() <<
") is less than ";
844 else MELAerr << nRequested_AssociatedPhotons;
854 for (
auto it=
daughters.cbegin()+1; it!=
daughters.cend(); it++){ firstPart.second = firstPart.second + it->second; }
864 costhetastar, costheta1, costheta2, Phi, Phi1, mVstar, mV,
869 aparts.at(0).second, aparts.at(0).first,
870 aparts.at(1).second, aparts.at(1).first,
876 if (!std::isfinite(costhetastar)) costhetastar=0;
877 if (!std::isfinite(costheta1)) costheta1=0;
878 if (!std::isfinite(costheta2)) costheta2=0;
879 if (!std::isfinite(Phi)) Phi=0;
880 if (!std::isfinite(Phi1)) Phi1=0;
881 if (!std::isfinite(mVstar)) mVstar=0;
882 if (!std::isfinite(mV)) mV=0;
885 <<
"Mela::computeVHAngles: (mVstar, mV, h1, h2, Phi, hs, Phi1) = "
886 << mVstar <<
", " << mV <<
", " << costheta1 <<
", " << costheta2 <<
", " << Phi <<
", "
887 << costhetastar <<
", " << Phi1 << endl;
928 =costheta1=costheta2=Phi=costhetastar=Phi1
929 =hbb=hWW=Phibb=Phi1bb
931 =hWminusf=PhiWminusf=0;
935 TLorentzVector
const nullVector(0, 0, 0, 0);
938 int nRequested_Tops=1;
939 int nRequested_Antitops=1;
943 simple_event_record mela_event;
944 mela_event.AssociationCode=partIncCode;
945 mela_event.nRequested_Tops=nRequested_Tops;
946 mela_event.nRequested_Antitops=nRequested_Antitops;
957 for (
unsigned int itd=0; itd<mela_event.pTopDaughters.at(0).size(); itd++) topDaughters.push_back(mela_event.pTopDaughters.at(0).at(itd));
958 for (
unsigned int itd=0; itd<mela_event.pAntitopDaughters.at(0).size(); itd++) antitopDaughters.push_back(mela_event.pAntitopDaughters.at(0).at(itd));
961 for (
unsigned int itop=0; itop<mela_event.pStableTops.size(); itop++) topDaughters.push_back(mela_event.pStableTops.at(itop));
962 for (
unsigned int itop=0; itop<mela_event.pStableAntitops.size(); itop++) antitopDaughters.push_back(mela_event.pStableAntitops.at(itop));
965 if (topDecay==0 && (mela_event.pStableTops.size()<1 || mela_event.pStableAntitops.size()<1)){
967 <<
"TUtil::TTHiggsMatEl: Number of stable tops (" << mela_event.pStableTops.size() <<
")"
968 <<
" and number of stable antitops (" << mela_event.pStableAntitops.size() <<
")"
969 <<
" in ttH process are not 1!" << endl;
972 else if (topDecay>0 && (mela_event.pTopDaughters.size()<1 || mela_event.pAntitopDaughters.size()<1)){
974 <<
"Mela::computeTTHAngles: Number of set of top daughters (" << mela_event.pTopDaughters.size() <<
")"
975 <<
" and number of set of antitop daughters (" << mela_event.pAntitopDaughters.size() <<
")"
976 <<
" in ttH process are not 1!" << endl;
979 else if (topDecay>0 && (mela_event.pTopDaughters.at(0).size()!=3 || mela_event.pAntitopDaughters.at(0).size()!=3)){
981 <<
"Mela::computeTTHAngles: Number of top daughters (" << mela_event.pTopDaughters.at(0).size() <<
")"
982 <<
" and number of antitop daughters (" << mela_event.pAntitopDaughters.at(0).size() <<
")"
983 <<
" in ttH process are not 3!" << endl;
986 if (topDaughters.size()<3){
for (
size_t ip=topDaughters.size(); ip<3; ip++) topDaughters.emplace_back(-9000, nullVector); }
987 if (antitopDaughters.size()<3){
for (
size_t ip=antitopDaughters.size(); ip<3; ip++) antitopDaughters.emplace_back(-9000, nullVector); }
993 for (
auto it=
daughters.cbegin()+1; it!=
daughters.cend(); it++){ firstPart.second = firstPart.second + it->second; }
1006 for (
size_t ip=0; ip<topDaughters.size(); ip++){
1007 if (ip>0){ pT += topDaughters.at(ip).second; pW += topDaughters.at(ip).second; }
1008 else pT += topDaughters.at(ip).second;
1016 for (
size_t ip=0; ip<antitopDaughters.size(); ip++){
1017 if (ip>0){ pT += antitopDaughters.at(ip).second; pW += antitopDaughters.at(ip).second; }
1018 else pT += antitopDaughters.at(ip).second;
1025 costhetastar, costheta1, costheta2, Phi, Phi1,
1026 hbb, hWW, Phibb, Phi1bb,
1028 hWminusf, PhiWminusf,
1035 topDaughters.at(0).second, topDaughters.at(0).first,
1036 topDaughters.at(1).second, topDaughters.at(1).first,
1037 topDaughters.at(2).second, topDaughters.at(2).first,
1039 antitopDaughters.at(0).second, antitopDaughters.at(0).first,
1040 antitopDaughters.at(1).second, antitopDaughters.at(1).first,
1041 antitopDaughters.at(2).second, antitopDaughters.at(2).first,
1048 if (!std::isfinite(costhetastar)) costhetastar=0;
1049 if (!std::isfinite(costheta1)) costheta1=0;
1050 if (!std::isfinite(costheta2)) costheta2=0;
1051 if (!std::isfinite(Phi)) Phi=0;
1052 if (!std::isfinite(Phi1)) Phi1=0;
1054 if (!std::isfinite(hbb)) hbb=0;
1055 if (!std::isfinite(hWW)) hWW=0;
1056 if (!std::isfinite(Phibb)) Phibb=0;
1057 if (!std::isfinite(Phi1bb)) Phi1bb=0;
1059 if (!std::isfinite(hWplusf)) hWplusf=0;
1060 if (!std::isfinite(PhiWplusf)) PhiWplusf=0;
1062 if (!std::isfinite(hWminusf)) hWminusf=0;
1063 if (!std::isfinite(PhiWminusf)) PhiWminusf=0;
1066 <<
"Mela::computeTTHAngles: (h1, h2, Phi, hs, Phi1) = "
1067 << costheta1 <<
", " << costheta2 <<
", " << Phi <<
", "
1068 << costhetastar <<
", " << Phi1 << endl;
1099 for (
int im=0; im<2; im++){
1101 selfDHzzcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1102 selfDHwwcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1113 double selfDZqqcoupl_input[
SIZE_ZQQ][2],
1114 double selfDZvvcoupl_input[
SIZE_ZVV][2],
1118 for (
int im=0; im<2; im++){
1129 double selfDZvvcoupl_input[
SIZE_ZVV][2],
1137 for (
int im=0; im<2; im++){
1147 double selfDGggcoupl_input[
SIZE_GGG][2],
1148 double selfDGqqcoupl_input[
SIZE_GQQ][2],
1149 double selfDGvvcoupl_input[
SIZE_GVV][2],
1153 for (
int im=0; im<2; im++){
1165 double selfDGggcoupl_input[
SIZE_GGG][2],
1166 double selfDGvvcoupl_input[
SIZE_GVV][2],
1173 for (
int im=0; im<2; im++){
1192 TLorentzVector
const nullVector(0, 0, 0, 0);
1193 float mZZ=0, mZ1=0, mZ2=0, costheta1=0, costheta2=0, Phi=0, costhetastar=0, Phi1=0;
1198 costheta1, costheta2, Phi,
1209 Y_rrv->setConstant(
true);
1212 if (computeAnaMELA){
1215 prob = integral->getVal();
1218 else prob =
pdf->getVal();
1222 Y_rrv->setConstant(
false);
1235 costheta1, costheta2, Phi,
1240 MELAout <<
"Mela::computeP: Condition (myME_ == TVar::MCFM && myProduction_ == TVar::ZZINDEPENDENT && myModel_ == TVar::bkgZZ/WW/ZGamma/ZJJ/GammaGamma)." << endl;
1241 vector<TLorentzVector> pDauVec =
calculate4Momentum(mZZ, mZ1, mZ1, acos(costhetastar), acos(costheta1), acos(costheta2), Phi1, Phi);
1243 <<
"\tOriginal mZZ=" << mZZ <<
" "
1244 <<
"m1=" << mZ1 <<
" "
1245 <<
"m2=" << mZ2 <<
" "
1246 <<
"h1=" << costheta1 <<
" "
1247 <<
"h2=" << costheta2 <<
" "
1248 <<
"Phi=" << Phi <<
" "
1249 <<
"hs=" << costhetastar <<
" "
1250 <<
"Phi1=" << Phi1 << endl;
1251 MELAout <<
"\tfor daughters:" << endl;
1252 for (
int iv=0; iv<2; iv++){
1256 <<
"x=" << pDauVec.at(2*iv+idau).X() <<
" "
1257 <<
"y=" << pDauVec.at(2*iv+idau).Y() <<
" "
1258 <<
"z=" << pDauVec.at(2*iv+idau).Z() <<
" "
1259 <<
"t=" << pDauVec.at(2*iv+idau).T() << endl;
1266 int gridsize_hs = 5;
1269 double hs_step = (hs_max - hs_min) /
double(gridsize_hs);
1271 int gridsize_phi1 = 5;
1272 double phi1_min = 0;
1273 double phi1_max = TMath::Pi();
1274 double phi1_step = (phi1_max - phi1_min) /
double(gridsize_phi1);
1276 for (
int i_hs = 0; i_hs < gridsize_hs + 1; i_hs++){
1277 double hs_val = hs_min + i_hs * hs_step;
1278 for (
int i_phi1 = 0; i_phi1 < gridsize_phi1 + 1; i_phi1++){
1279 double phi1_val = phi1_min + i_phi1 * phi1_step;
1284 vector<TLorentzVector> pDauVec =
calculate4Momentum(mZZ, mZ1, mZ2, acos(hs_val), acos(costheta1), acos(costheta2), phi1_val, Phi);
1285 for (
int iv=0; iv<2; iv++){
1292 MELAout <<
"Mela::computeP: hs, Phi1 are now " << hs_val <<
" " << phi1_val << endl;
1293 unsigned int idau=1;
1295 MELAout <<
"Dau " << idau <<
" "
1296 <<
"id=" << tmpPart.first <<
" "
1297 <<
"x=" << tmpPart.second.X() <<
" "
1298 <<
"y=" << tmpPart.second.Y() <<
" "
1299 <<
"z=" << tmpPart.second.Z() <<
" "
1300 <<
"t=" << tmpPart.second.T() << endl;
1304 vector<MELAParticle*> partList_tmp;
1305 vector<MELACandidate*> candList_tmp;
1321 for (
MELAParticle* tmpPart:partList_tmp)
delete tmpPart;
1326 prob = prob / float((gridsize_hs + 1) * (gridsize_phi1 +1));
1351 coupl_mix[0][0][0] =1.;
1352 coupl_mix[0][3][0] =2.521;
1353 coupl_1[0][0][0] =1.;
1354 coupl_2[0][3][0] =2.521;
1357 coupl_mix[0][0][0] =1.;
1358 coupl_mix[0][3][1] =2.521;
1359 coupl_1[0][0][0] =1.;
1360 coupl_2[0][3][1] =2.521;
1363 coupl_mix[0][0][0] =1.;
1364 coupl_mix[0][1][0] = 1.638;
1365 coupl_1[0][0][0] =1.;
1366 coupl_2[0][1][0] = 1.638;
1369 coupl_mix[0][0][0] =1.;
1370 coupl_mix[0][1][1] = 1.638;
1371 coupl_1[0][0][0] =1.;
1372 coupl_2[0][1][1] = 1.638;
1375 coupl_mix[0][0][0] =1.;
1376 coupl_mix[0][11][0] = 12046.01;
1377 coupl_1[0][0][0] =1.;
1378 coupl_2[0][11][0] = 12046.01;
1381 coupl_mix[0][0][0] =1.;
1382 coupl_mix[0][4][0] = 0.0688;
1383 coupl_1[0][0][0] =1.;
1384 coupl_2[0][4][0] = 0.0688;
1387 coupl_mix[0][0][0] =1.;
1388 coupl_mix[0][7][0] = -0.0898;
1389 coupl_1[0][0][0] =1.;
1390 coupl_2[0][7][0] = -0.0898;
1393 coupl_mix[0][0][0] =1.;
1394 coupl_mix[0][6][0] = 0.0855;
1395 coupl_1[0][0][0] =1.;
1396 coupl_2[0][6][0] = 0.0855;
1399 coupl_mix[0][0][0] =1.;
1400 coupl_mix[0][9][0] = -0.0907;
1401 coupl_1[0][0][0] =1.;
1402 coupl_2[0][9][0] = -0.0907;
1405 coupl_mix[0][0][0] =1.;
1406 coupl_mix[0][30][0] = -7591.914;
1407 coupl_1[0][0][0] =1.;
1408 coupl_2[0][30][0] = -7591.914;
1411 coupl_mix[0][0][0] =1.;
1412 coupl_mix[0][30][1] = -7591.914;
1413 coupl_1[0][0][0] =1.;
1414 coupl_2[0][30][1] = -7591.914;
1417 MELAout <<
"Error: Not supported!"<<endl;
1425 prob = pMix- p1- p2;
1434 double selfDAZffcoupl_input[
SIZE_AZff][2],
1439 for (
int im=0; im<2; im++){
1441 selfDHzzcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1442 selfDHwwcoupl[jh][ic][im] = selfDHwwcoupl_input[jh][ic][im];
1446 for (
int im=0; im<2; im++){
1449 for (
int im=0; im<2; im++){
1465 bool hasFailed =
false;
1467 MELAout <<
"Mela::computeProdDecP ME is not supported for ME " <<
myME_ << endl;
1483 MELAout <<
"Mela::computeProdDecP production mode is not supported for production " <<
myProduction_ << endl;
1487 if (hasFailed) prob=0;
1502 double selfDHggcoupl_input[
SIZE_HGG][2],
1508 for (
int im=0; im<2; im++){
for (
int ic=0; ic<
SIZE_HGG; ic++)
selfDHggcoupl[0][ic][im] = selfDHggcoupl_input[ic][im]; }
1510 for (
int im=0; im<2; im++){
1512 selfDHzzcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1513 selfDHwwcoupl[jh][ic][im] = selfDHwwcoupl_input[jh][ic][im];
1534 TLorentzVector nullFourVector(0, 0, 0, 0);
1535 bool isJet2Fake =
false;
1538 unsigned int firstJetIndex=0;
1539 TLorentzVector jet1, higgs;
1540 TLorentzVector jet1massless(0, 0, 0, 0);
1541 TLorentzVector jet2massless(0, 0, 0, 0);
1548 if (tmpPart->passSelection){
1564 const double threshold = 1000.*
LHCsqrts/2.;
1565 TLorentzVector pTotal = higgs+jet1massless+jet2massless;
1566 double sysZ = pTotal.Z();
1567 if (fabs(sysZ)>threshold){
1568 double maxpz2 = threshold - higgs.Z() - jet1massless.Z();
1569 if (fabs(maxpz2)>0.){
1570 double ratio = jet2massless.Z()/maxpz2;
1571 double absp=sqrt(pow(jet2massless.Pt(), 2)+pow(jet2massless.Z()*
ratio, 2));
1572 if (
myVerbosity_>=
TVar::INFO)
MELAout <<
"Mela::computeProdP, isJet2Fake=true case: Rescaling pz of fake jet by " <<
ratio <<
" and energy = " << absp <<
"." << endl;
1573 jet2massless.SetXYZT(jet2massless.X(), jet2massless.Y(), jet2massless.Z()*
ratio, absp);
1576 if (
myVerbosity_>=
TVar::INFO)
MELAout <<
"Mela::computeProdP, isJet2Fake=true case: Unable to rescaling pz of fake jet since max(|pz|)<0. Setting to 0 with appropriate energy = pT = " << jet2massless.Pt() <<
"." << endl;
1577 jet2massless.SetXYZT(jet2massless.X(), jet2massless.Y(), 0., jet2massless.Pt());
1586 firstJet->
p4.SetXYZT(jet1massless.X(), jet1massless.Y(), jet1massless.Z(), jet1massless.T());
1595 std::vector<double> etaArray;
1596 std::vector<double> pArray;
1597 double eta_max = 10;
1598 if (jet2massless.Pt()>0.) eta_max = max(eta_max, 1.2*fabs(jet2massless.Eta()));
1599 double eta_min = -eta_max;
1601 for (
int iter=0; iter<nGrid; iter++){
1604 double jet2temp_eta = ((double)iter)*(eta_max-eta_min) / (((
double)nGrid) - 1.) + eta_min;
1605 etaArray.push_back(jet2temp_eta);
1606 double jet2temp_sinh_eta = TMath::SinH(jet2temp_eta);
1607 double jet2temp_pz = jet2massless.Pt()*jet2temp_sinh_eta;
1608 fakeJet.
p4.SetZ(jet2temp_pz);
1609 fakeJet.
p4.SetX(jet2massless.X()); fakeJet.
p4.SetY(jet2massless.Y()); fakeJet.
p4.SetT(fakeJet.
p4.P());
1612 const double threshold = 1000.*
LHCsqrts/2.;
1613 TLorentzVector pTotal = higgs+jet1massless+fakeJet.
p4;
1614 double sys = (pTotal.T()+fabs(pTotal.Z()))/2.;
1615 if (fabs(sys)<threshold){
1619 pArray.push_back((
double)prob_temp);
1622 const double grid_precision = 0.15;
1624 for (
int iG=0; iG<nGrid-1; iG++){
1625 if (pArray[iG]==pArray[iG+1])
continue;
1626 if (etaArray[iG]==etaArray[iG+1])
continue;
1630 TGraph interpolator(nGrid, etaArray.data(), pArray.data());
1631 double derivative_first = (pArray[1]-pArray[0])/(etaArray[1]-etaArray[0]);
1632 double derivative_last = (pArray[nGrid-1]-pArray[nGrid-2])/(etaArray[nGrid-1]-etaArray[nGrid-2]);
1633 TSpline3 spline(
"spline", &interpolator,
"b1e1", derivative_first, derivative_last);
1634 double x_middle = (etaArray[iG]+etaArray[iG+1])*0.5;
1635 double y_middle = (pArray[iG]+pArray[iG+1])*0.5;
1636 double y_sp = spline.Eval(x_middle);
1637 if (y_sp<0) y_sp = 0;
1639 std::vector<double>::iterator gridIt;
1641 if (fabs(y_sp-y_middle)<grid_precision*fabs(y_middle) || fabs(etaArray[iG+1]-etaArray[iG])<1e-3){
1642 gridIt = pArray.begin()+iG+1;
1643 pArray.insert(gridIt, y_sp);
1644 gridIt = etaArray.begin()+iG+1;
1645 etaArray.insert(gridIt, x_middle);
1651 double jet2temp_eta = x_middle;
1652 gridIt = etaArray.begin()+iG+1;
1653 etaArray.insert(gridIt, x_middle);
1654 double jet2temp_sinh_eta = TMath::SinH(jet2temp_eta);
1655 double jet2temp_pz = jet2massless.Pt()*jet2temp_sinh_eta;
1656 fakeJet.
p4.SetZ(jet2temp_pz);
1657 fakeJet.
p4.SetX(jet2massless.X()); fakeJet.
p4.SetY(jet2massless.Y()); fakeJet.
p4.SetT(fakeJet.
p4.P());
1660 const double threshold = 1000.*
LHCsqrts/2.;
1661 TLorentzVector pTotal = higgs+jet1massless+fakeJet.
p4;
1662 double sys = (pTotal.T()+fabs(pTotal.Z()))/2.;
1663 if (fabs(sys)<threshold){
1667 gridIt = pArray.begin()+iG+1;
1668 pArray.insert(gridIt, (
double)prob_temp);
1677 int iGFirst=0, iGLast=nGrid-1;
1678 for (
int iG=1; iG<nGrid; iG++){
1679 if (pArray[iG]>0 && pArray[iG-1]==0){
1684 for (
int iG=nGrid-2; iG>=0; iG--){
1685 if (pArray[iG]>0 && pArray[iG+1]==0){
1690 double dEtaGrid = etaArray[iGLast] - etaArray[iGFirst];
1691 for (
int iG=iGFirst; iG<iGLast-1; iG++){
1692 double dEta = etaArray[iG+1] - etaArray[iG];
1693 double sumProb = pArray[iG]+pArray[iG+1];
1695 dEta = dEta/dEtaGrid;
1696 double addProb = sumProb*dEta;
1700 firstJet->
p4.SetXYZT(jet1.X(), jet1.Y(), jet1.Z(), jet1.T());
1705 if (
melaCand!=candOriginal)
MELAerr <<
"Mela::computeProdP: melaCand!=candOriginal at the end of the fake jet scenario!" << endl;
1733 bool includeHiggsDecay,
1739 for (
int im=0; im<2; im++){
1742 selfDHzzcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1743 selfDHwwcoupl[jh][ic][im] = selfDHvvcoupl_input[jh][ic][im];
1756 bool includeHiggsDecay,
1812 bool hasFailed=
false;
1813 int id_original[2][2];
1814 for (
int iv=0; iv<2; iv++){
1818 for (
int ivd=0; ivd<2; ivd++) id_original[iv][ivd]=Vi->
getDaughter(ivd)->
id;
1845 float dXsec_HZZ_JHU, dXsec_HZZ_JHU_interf;
1850 computeP(dXsec_HZZ_JHU_interf,
false);
1853 w = dXsec_HZZ_JHU_interf / dXsec_HZZ_JHU;
1870 bool hasFailed=
false;
1871 int id_original[2][2];
1872 for (
int iv=0; iv<2; iv++){
1876 for (
int ivd=0; ivd<2; ivd++) id_original[iv][ivd]=Vi->
getDaughter(ivd)->
id;
1880 if (abs(id_original[0][0])==11 && abs(id_original[1][0])==11 && abs(id_original[0][1])==11 && abs(id_original[1][1])==11)
super->
SetDecayChannel(
"4e");
1881 else if (abs(id_original[0][0])==13 && abs(id_original[1][0])==13 && abs(id_original[0][1])==13 && abs(id_original[1][1])==13)
super->
SetDecayChannel(
"4mu");
1883 (abs(id_original[0][0])==11 && abs(id_original[0][1])==11 && abs(id_original[1][0])==13 && abs(id_original[1][1])==13)
1885 (abs(id_original[0][0])==13 && abs(id_original[0][1])==13 && abs(id_original[1][0])==11 && abs(id_original[1][1])==11)
1887 else{
if (
myVerbosity_>=
TVar::ERROR)
MELAerr <<
"Mela::computePM4l: SuperMELA is currently not implemented for decay states other than 4e. 4mu, 2e2mu." << endl; hasFailed=
true; }
1908 std::pair<double, double> m4lP =
super->
M4lProb(mZZtmp);
1921 float bkg_VAMCFM_noscale,
1922 float ggzz_VAMCFM_noscale,
1923 float ggHZZ_prob_pure_noscale,
1924 float ggHZZ_prob_int_noscale,
1928 float total_sig_ME = (widthScale * ggHZZ_prob_pure_noscale + sqrt(widthScale) * ggHZZ_prob_int_noscale + ggzz_VAMCFM_noscale);
1929 float total_bkg_ME = bkg_VAMCFM_noscale;
1930 float kd_denominator = (total_sig_ME+total_bkg_ME);
1931 if (kd_denominator>0.) myDggr = total_sig_ME/(total_sig_ME+total_bkg_ME);
1941 MELAout <<
"Only support MCFM and D_gg10"<<endl;
1947 float bkg_VAMCFM, ggzz_VAMCFM_noscale, ggHZZ_prob_pure_noscale, ggHZZ_prob_int_noscale, bkgHZZ_prob_noscale;
1953 bkgHZZ_prob_noscale /= ggScale;
1954 ggHZZ_prob_pure_noscale /= ggScale;
1955 ggzz_VAMCFM_noscale /= ggScale;
1957 ggHZZ_prob_int_noscale = bkgHZZ_prob_noscale - ggHZZ_prob_pure_noscale - ggzz_VAMCFM_noscale;
1961 constructDggr(bkg_VAMCFM, ggzz_VAMCFM_noscale, ggHZZ_prob_pure_noscale, ggHZZ_prob_int_noscale, 10., prob);
2012 for (
int im=0; im<2; im++){
2136 Double_t fppReal = 1./sqrt(6.) * (
c1/4.*2. + 2.*
c2);
2137 Double_t fppImag = 1./sqrt(6.) *
c5;
2138 Double_t fmmReal = 1./sqrt(6.) * (
c1/4.*2. + 2.*
c2);
2139 Double_t fmmImag = 1./sqrt(6.)*
c5;
2140 Double_t fmpReal = 1./4.*
c1*2.;
2141 Double_t fmpImag = 0;
2142 Double_t fpp = fppImag*fppImag + fppReal*fppReal;
2143 Double_t fmm = fmmImag*fmmImag + fmmReal*fmmReal;
2144 Double_t fmp = fmpImag*fmpImag + fmpReal*fmpReal;
2157 MELAout <<
"Mela::configureAnalyticalPDFs -> ERROR TVar::Process not applicable!!! ME: " <<
myME_ <<
", model: " <<
myModel_ << endl;
2161 if (
pdf==0) noPass=
true;
2246 if (
std::isnan(constant) || std::isinf(constant) || constant<=0.) constant=0;
2247 else constant=1./constant;
2256 unsigned int iarray=0;
2281 else return constant;
2309 const int decid = 121;
2316 const bool is4mu = (decid==28561);
2317 const bool is4e = (decid==14641 || decid==50625);
2318 const bool is2mu2e = (decid==20449 || decid==27225 || decid==38025 || decid==169 || decid==121);
2320 const unsigned int nPossibleHandles=6;
2323 float constant_tmp=0;
2433 if (hasVBF) pchandle[0]=hvbf;
2434 if (hasZH) pchandle[1]=hzh;
2435 if (hasWH) pchandle[2]=hwh;
2438 if (hasVBF) pchandle[3]=hvbs;
2439 if (hasZH) pchandle[4]=hzzz;
2440 if (hasWH) pchandle[5]=hwzz;
2461 bool hasNullHandle=
true;
2462 for (
unsigned int ihandle=0; ihandle<nPossibleHandles; ihandle++){
if (pchandle[ihandle]!=0){ constant_tmp += pchandle[ihandle]->
Eval(
getIORecord(),
myVerbosity_); hasNullHandle=
false; } }
2463 if (hasNullHandle)
return constant;
2465 constant = constant_tmp;
2529 TString filename, spname;
2532 filename =
"pAvgSmooth_JHUGen_JJQCD_HSMHiggs";
2533 spname =
"P_ConserveDifermionMass";
2535 spname =
"P_MomentumToEnergy";
2539 filename =
"pAvgSmooth_JHUGen_JQCD_HSMHiggs";
2540 spname =
"P_ConserveDifermionMass";
2542 spname =
"P_MomentumToEnergy";
2546 filename =
"pAvgSmooth_JHUGen_JJVBF_HSMHiggs";
2547 spname =
"P_ConserveDifermionMass";
2549 spname =
"P_MomentumToEnergy";
2553 filename =
"pAvgSmooth_JHUGen_Had_ZH_HSMHiggs";
2554 spname =
"P_ConserveDifermionMass";
2556 spname =
"P_MomentumToEnergy";
2560 filename =
"pAvgSmooth_JHUGen_Had_WH_HSMHiggs";
2561 spname =
"P_ConserveDifermionMass";
2563 spname =
"P_MomentumToEnergy";
2568 filename =
"pAvgSmooth_MCFM_JJQCD_bkgZJets_13TeV_2l2q";
2569 spname =
"P_ConserveDifermionMass";
2573 filename =
"pAvgSmooth_JHUGen_ZZGG_HSMHiggs";
2574 spname =
"P_ConserveDifermionMass_4mu";
2576 spname =
"P_ConserveDifermionMass_4e";
2578 spname =
"P_ConserveDifermionMass_2mu2e";
2581 filename =
"pAvgSmooth_MCFM_ZZGG_HSMHiggs";
2582 spname =
"P_ConserveDifermionMass_4mu";
2584 spname =
"P_ConserveDifermionMass_4e";
2586 spname =
"P_ConserveDifermionMass_2mu2e";
2589 filename =
"pAvgSmooth_MCFM_JJVBF_S_HSMHiggs";
2590 spname =
"P_ConserveDifermionMass_4mu";
2592 spname =
"P_ConserveDifermionMass_4e";
2594 spname =
"P_ConserveDifermionMass_2mu2e";
2597 filename =
"pAvgSmooth_MCFM_Had_ZH_S_HSMHiggs";
2598 spname =
"P_ConserveDifermionMass_4mu";
2600 spname =
"P_ConserveDifermionMass_4e";
2602 spname =
"P_ConserveDifermionMass_2mu2e";
2605 filename =
"pAvgSmooth_MCFM_Had_WH_S_HSMHiggs";
2606 spname =
"P_ConserveDifermionMass_4mu";
2608 spname =
"P_ConserveDifermionMass_4e";
2610 spname =
"P_ConserveDifermionMass_2mu2e";
2614 filename =
"pAvgSmooth_MCFM_ZZGG_bkgZZ";
2615 spname =
"P_ConserveDifermionMass_4mu";
2617 spname =
"P_ConserveDifermionMass_4e";
2619 spname =
"P_ConserveDifermionMass_2mu2e";
2622 filename =
"pAvgSmooth_MCFM_ZZQQB_bkgZZ";
2623 spname =
"P_ConserveDifermionMass_4mu";
2625 spname =
"P_ConserveDifermionMass_4e";
2627 spname =
"P_ConserveDifermionMass_2mu2e";
2630 filename =
"pAvgSmooth_MCFM_JJVBF_bkgZZ";
2631 spname =
"P_ConserveDifermionMass_4mu";
2633 spname =
"P_ConserveDifermionMass_4e";
2635 spname =
"P_ConserveDifermionMass_2mu2e";
2638 filename =
"pAvgSmooth_MCFM_Had_ZH_bkgZZ";
2639 spname =
"P_ConserveDifermionMass_4mu";
2641 spname =
"P_ConserveDifermionMass_4e";
2643 spname =
"P_ConserveDifermionMass_2mu2e";
2646 filename =
"pAvgSmooth_MCFM_Had_WH_bkgZZ";
2647 spname =
"P_ConserveDifermionMass_4mu";
2649 spname =
"P_ConserveDifermionMass_4e";
2651 spname =
"P_ConserveDifermionMass_2mu2e";
2654 filename =
"pAvgSmooth_MCFM_JJQCD_bkgZZ";
2655 spname =
"P_ConserveDifermionMass_4mu";
2657 spname =
"P_ConserveDifermionMass_4e";
2659 spname =
"P_ConserveDifermionMass_2mu2e";
2676 string cfile_fullpath;
2681 const string path = MELAPKGPATH +
"data/";
2685 const unsigned int npossiblesqrts=3;
2686 const double possible_sqrts[npossiblesqrts]={ 7, 8, 13 };
2687 vector<double> trysqrts;
2688 for (
unsigned isq=0; isq<npossiblesqrts; isq++){
2689 double val = possible_sqrts[isq];
2691 bool inserted=
false;
2692 for (std::vector<double>::iterator it = trysqrts.begin(); it<trysqrts.end(); it++){
2695 trysqrts.insert(it, val);
2699 if (!inserted) trysqrts.push_back(val);
2701 for (
auto& dsqrts:trysqrts){
2702 TString strsqrts = Form(
"%s_%.0f%s", relpath.Data(), dsqrts,
"TeV");
2703 cfile_fullpath = path;
2704 cfile_fullpath.append(strsqrts.Data());
2705 cfile_fullpath.append(
".root");
2706 pchandle =
new MelaPConstant(me_, prod_, proc_, cfile_fullpath.c_str(), spname.Data());
2708 if (
myVerbosity_>=
TVar::DEBUG)
MELAout <<
"Mela::getPConstantHandle: Full path and spline name: " << cfile_fullpath <<
", " << spname <<
" is valid." << endl;
2712 if (
myVerbosity_>=
TVar::DEBUG)
MELAout <<
"Mela::getPConstantHandle: Full path and spline name: " << cfile_fullpath <<
", " << spname <<
" is invalid." << endl;
2718 cfile_fullpath = path;
2719 cfile_fullpath.append(relpath.Data());
2720 cfile_fullpath.append(
".root");
2721 pchandle =
new MelaPConstant(me_, prod_, proc_, cfile_fullpath.c_str(), spname.Data());
2725 if (
myVerbosity_>=
TVar::DEBUG && pchandle==0)
MELAerr <<
"Mela::getPConstantHandle: Handle of " << spname <<
" from " << cfile_fullpath <<
" is invalid!" << endl;
2787 delete handle; handle=0;
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_S_HSMHiggs_4e
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_bkgZZ_4e
RooAbsReal * gzgs4List[1][2]
void reset_Width(double inmass, int ipart)
void getConstant(float &prob)
This function returns a multiplicative constant to the matrix element calculation.
std::pair< int, TLorentzVector > SimpleParticle_t
void scaleMomentumToEnergy(const TLorentzVector &massiveJet, TLorentzVector &masslessJet, double mass=0)
RooSpinTwo::modelCouplings couplings
void computeP_selfDspin2(double selfDGggcoupl_input[SIZE_GGG][2], double selfDGqqcoupl_input[SIZE_GQQ][2], double selfDGvvcoupl_input[SIZE_GVV][2], float &prob, bool useConstant=true)
This is a function that calls Mela::computeP with preset quark, gluon, and vector boson couplings for...
bool isALepton(const int id)
float getConstant_FourFermionDecay(const int &decid)
MelaPConstant * pAvgSmooth_MCFM_ZZQQB_bkgZZ_4mu
void setVerbosity(TVar::VerbosityLevel verbosity)
MELAParticle * getDaughter(int index) const
MelaPConstant * pAvgSmooth_JHUGen_JQCD_HSMHiggs[TVar::nFermionMassRemovalSchemes-1]
MelaPConstant * pAvgSmooth_MCFM_JJQCD_bkgZZ_4e
void setInputEvent(SimpleParticleCollection_t *pDaughters, SimpleParticleCollection_t *pAssociated=0, SimpleParticleCollection_t *pMothers=0, bool isGen=false)
Sets the input event for MELA. MELA cannot run without this.
virtual void resetHypotheses()
double selfDGqqcoupl[SIZE_GQQ][2]
double selfDGvpvpcoupl[SIZE_GVV][2]
void getPConstantHandles()
virtual void setTensorPolarization(int ig, double initval)
void build(double mh_)
This is the actual building of the tool that occurs in each instance of the Mela::Mela constructor.
void addAssociatedJet(MELAParticle *myParticle)
void get_XPropagator(TVar::ResonancePropagatorScheme scheme, float &prop)
void resetMCFM_EWKParameters(double ext_Gf, double ext_aemmz, double ext_mW, double ext_mZ, double ext_xW, int ext_ewscheme=3)
Resets the electroweak parameters back to their defaults.
double getPrimaryWidth(int ipart)
A function to get the current primary EW/QCD parameters from MELA.
int getNCandidates()
Returns the size of the candidate list TEvtProb::candList.
void computeDecayAngles(float &qH, float &m1, float &m2, float &costheta1, float &costheta2, float &Phi, float &costhetastar, float &Phi1)
computes the decay angles for gg -> H -> ZZ as defined in Figure 1 of arXiv:1001.3396
MelaPConstant * getPConstantHandle(TVar::MatrixElement me_, TVar::Production prod_, TVar::Process proc_, TString relpath, TString spname, const bool useSqrts=false)
void computeProdP_VH(double selfDHvvcoupl_input[nSupportedHiggses][SIZE_HVV][2], float &prob, bool includeHiggsDecay=false, bool useConstant=true)
void set_mHiggs_wHiggs(double mh_, double gah_, int index)
virtual void addHypothesis(int ig, double initval, double iphase=0)
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_bkgZZ_2mu2e
void compute4FermionWeight(float &w)
int getNDaughters() const
MelaPConstant * pAvgSmooth_MCFM_ZZGG_HSMHiggs_4e
std::vector< MELATopCandidate_t * > * get_TopCandidateCollection()
int selfDHwwCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ]
void reset_Mass(double inmass, int ipart)
MelaPConstant * pAvgSmooth_MCFM_JJVBF_S_HSMHiggs_4mu
@ nFermionMassRemovalSchemes
MelaIO * getIORecord()
Returns the MELAIO object, and by consequence, the entire parton-by-parton matrix element record.
MelaPConstant * pAvgSmooth_MCFM_Had_WH_bkgZZ_4mu
MelaPConstant * pAvgSmooth_MCFM_Had_WH_S_HSMHiggs_2mu2e
@ kUseAssociated_UnstableTops
void setCurrentCandidateFromIndex(unsigned int icand)
Switches the candidate that you are working on to another candidate based off of an index.
MELAParticle * getSortedV(int index) const
TString ProductionName(TVar::Production temp)
MelaPConstant * pAvgSmooth_MCFM_ZZGG_bkgZZ_4mu
double selfDHttcoupl[nSupportedHiggses][SIZE_HQQ][2]
MelaPConstant * pAvgSmooth_MCFM_JJQCD_bkgZZ_2mu2e
void setRemoveLeptonMasses(bool MasslessLeptonSwitch=true)
either permits or forbids massive leptons.
void getPAux(float &prob)
void SetVerbosity(bool verb=true)
RooSpin::modelMeasurables getMeasurablesRRV()
Returns a RooSpin::modelMeasureables object containing many observable quantities.
void setVerbosity(TVar::VerbosityLevel verbosity_=TVar::ERROR)
Sets the verbosity for MELA outside of the initial constructor.
void setCurrentCandidate(MELACandidate *cand)
Switches the candidate that you are working on to another candidate object specified.
RooAbsReal * g2List[8][2]
RooAbsReal * gzgs3List[1][2]
void resetWidth(double inwidth, int ipart)
Resets the width for a particle that is an electroweak parameter according to its id.
void setMelaHiggsWidth(double myHiggsWidth=-1, int index=0)
Sets the width of your chosen Higgs.
void setSpinZeroCouplings()
void computeProdXS_JH(float &mevalue)
std::vector< MELATopCandidate_t * > * getTopCandidateCollection()
Same as getNCandidates but specifically for Top Quark Candidates.
void computeD_CP(TVar::MatrixElement myME, TVar::Process myType, float &prob)
computes the value of D_CP
void computeTTHAngles(float &hs, float &hincoming, float &hTT, float &PhiTT, float &Phi1, float &hbb, float &hWW, float &Phibb, float &Phi1bb, float &hWplusf, float &PhiWplusf, float &hWminusf, float &PhiWminusf, TLorentzVector p4M11, int Z1_lept1Id, TLorentzVector p4M12, int Z1_lept2Id, TLorentzVector p4M21, int Z2_lept1Id, TLorentzVector p4M22, int Z2_lept2Id, TLorentzVector b, int bId, TLorentzVector Wplusf, int WplusfId, TLorentzVector Wplusfb, int WplusfbId, TLorentzVector bbar, int bbarId, TLorentzVector Wminusf, int WminusfId, TLorentzVector Wminusfb, int WminusfbId, TLorentzVector *injet1=0, int injet1Id=0, TLorentzVector *injet2=0, int injet2Id=0)
MelaPConstant * pAvgSmooth_MCFM_Had_WH_S_HSMHiggs_4mu
MelaPConstant * pAvgSmooth_JHUGen_Had_WH_HSMHiggs[TVar::nFermionMassRemovalSchemes-1]
void computeProdP(double selfDHggcoupl_input[SIZE_HGG][2], double selfDHvvcoupl_input[nSupportedHiggses][SIZE_HVV][2], double selfDHwwcoupl_input[nSupportedHiggses][SIZE_HVV][2], float &prob, bool useConstant=true)
computes Production side probabilities while taking in coupling arrays
MelaPConstant * pAvgSmooth_MCFM_ZZQQB_bkgZZ_4e
double selfDHzpzpcoupl[SIZE_HVV][2]
void setRenFacScaleMode(TVar::EventScaleScheme renormalizationSch, TVar::EventScaleScheme factorizationSch, double ren_sf, double fac_sf)
Sets the renormalization and the factorization schemes.
void applyLeptonMassCorrection(bool flag=true)
double getMEConst() const
void SetVerbosity(const TVar::VerbosityLevel verbosity_)
double selfDWpffcoupl[SIZE_Vpff][2]
void reset_SelfDCouplings()
double selfDZpffcoupl[SIZE_Vpff][2]
double selfDaTQGCcoupl[SIZE_ATQGC][2]
float getConstant_JHUGenUndecayed()
void set_Verbosity(TVar::VerbosityLevel verbosity_)
void set_SpinZeroContact(double selfDHzzpcoupl[SIZE_HVV][2], double selfDHzpzpcoupl[SIZE_HVV][2], double selfDHwwpcoupl[SIZE_HVV][2], double selfDHwpwpcoupl[SIZE_HVV][2])
MelaPConstant * pAvgSmooth_JHUGen_ZZGG_HSMHiggs_2mu2e
void set_SpinOneCouplings(double selfDZqqcoupl[SIZE_ZQQ][2], double selfDZvvcoupl[SIZE_ZVV][2])
void computeVBFAngles(float &Q2V1, float &Q2V2, float &costheta1, float &costheta2, float &Phi, float &costhetastar, float &Phi1)
computes the decay angles for Vector Boson Fusion (VBF) -> H -> ZZ as defined in Figure 1 of arXiv:12...
void PrintCandidateSummary(MELACandidate *cand)
bool isAZBoson(const int id)
void setProcess(TVar::Process myModel, TVar::MatrixElement myME, TVar::Production myProduction)
Sets the process, matrix element, and production that MELA is to use for this event....
double get_HiggsWidthAtPoleMass(double mass)
logical function isnan(x)
Mela(double LHCsqrts_=13., double mh_=125., TVar::VerbosityLevel verbosity_=TVar::ERROR)
the MELA constructor
std::vector< TLorentzVector > Calculate4Momentum(double Mx, double M1, double M2, double theta, double theta1, double theta2, double Phi1, double Phi)
RooqqZZ_JHU_ZgammaZZ_fast * qqZZmodel
void resetQuarkMasses()
Resets the masses of each quark to their original values.
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_S_HSMHiggs_4mu
double selfDHqqcoupl[nSupportedHiggses][SIZE_HQQ][2]
void setTempCandidate(SimpleParticleCollection_t *pDaughters, SimpleParticleCollection_t *pAssociated=0, SimpleParticleCollection_t *pMothers=0, bool isGen=false)
Sets a temporary MELA candidate, by setting melaCand in Xcal2 to a temporary candidate without pushin...
void set_SpinZeroCouplings(double selfDHggcoupl[nSupportedHiggses][SIZE_HGG][2], double selfDHg4g4coupl[nSupportedHiggses][SIZE_HGG][2], double selfDHqqcoupl[nSupportedHiggses][SIZE_HQQ][2], double selfDHbbcoupl[nSupportedHiggses][SIZE_HQQ][2], double selfDHttcoupl[nSupportedHiggses][SIZE_HQQ][2], double selfDHb4b4coupl[nSupportedHiggses][SIZE_HQQ][2], double selfDHt4t4coupl[nSupportedHiggses][SIZE_HQQ][2], double selfDHzzcoupl[nSupportedHiggses][SIZE_HVV][2], double selfDHwwcoupl[nSupportedHiggses][SIZE_HVV][2], double selfDHzzLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ], double selfDHwwLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ], int selfDHzzCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ], int selfDHwwCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ], bool diffHWW=false)
void deletePConstantHandle(MelaPConstant *&handle)
ScalarPdfFactory_HVV * ggSpin0Model
RooRealVar * costheta1_rrv
MelaPConstant * pAvgSmooth_MCFM_ZZGG_bkgZZ_2mu2e
void deletePConstantHandles()
void append_TopCandidate(SimpleParticleCollection_t *TopDaughters)
double selfDGvvcoupl[SIZE_GVV][2]
MELACandidate * shallowCopy()
bool isAPhoton(const int id)
void set_TempCandidate(SimpleParticleCollection_t *pDaughters, SimpleParticleCollection_t *pAssociated=0, SimpleParticleCollection_t *pMothers=0, bool isGen=false)
void set_CandidateDecayMode(TVar::CandidateDecayMode mode)
double get_PrimaryMass(int ipart)
void set_AZffCouplings(double selfDAZffcoupl[SIZE_AZff][2])
void set_SpinTwoCouplings(double selfDGqqcoupl[SIZE_GQQ][2], double selfDGggcoupl[SIZE_GGG][2], double selfDGvvcoupl[SIZE_GVV][2])
void computeP_selfDspin0(double selfDHvvcoupl_input[nSupportedHiggses][SIZE_HVV][2], float &prob, bool useConstant=true)
This is a function that calls Mela::computeP with preset quark, gluon, and vector boson couplings for...
int get_CurrentCandidateIndex()
void set_Process(TVar::Process process_, TVar::MatrixElement me_, TVar::Production production_)
MELACandidate * get_CurrentCandidate()
MELAOutputStreamer MELAout
MelaPConstant * pAvgSmooth_MCFM_ZZQQB_bkgZZ_2mu2e
void set_mHiggs(double mh_, int index)
MelaPConstant * pAvgSmooth_MCFM_Had_WH_bkgZZ_2mu2e
bool differentiate_HWW_HZZ
void set_VprimeContactCouplings(double selfDZpffcoupl[SIZE_Vpff][2], double selfDWpffcoupl[SIZE_Vpff][2], double M_Zprime, double Ga_Zprime, double M_Wprime, double Ga_Wprime)
double getHiggsWidthAtPoleMass(double mass)
Returns the width of the Higgs at a given pole mass as a calculation.
void SetPathToCards(std::string dirToCards)
virtual void addHypothesis(int ig, int ilam, double iphase=0, double altparam_fracval=0)
MELACandidate * ConvertVectorFormat(SimpleParticleCollection_t *pDaughters, SimpleParticleCollection_t *pAssociated, SimpleParticleCollection_t *pMothers, bool isGen, std::vector< MELAParticle * > *particleList, std::vector< MELACandidate * > *candList)
void computeTTHAngles(int topDecay, float &mT1, float &mW1, float &mT2, float &mW2, float &costheta1, float &costheta2, float &Phi, float &costhetastar, float &Phi1, float &hbb, float &hWW, float &Phibb, float &Phi1bb, float &hWplusf, float &PhiWplusf, float &hWminusf, float &PhiWminusf)
MelaPConstant * pAvgSmooth_JHUGen_ZZGG_HSMHiggs_4e
static void cleanLinkedFiles()
TVar::MatrixElement myME_
std::pair< double, double > M4lProb(double m4l)
void computeVHAngles(float &mVstar, float &mV, float &costheta1, float &costheta2, float &Phi, float &costhetastar, float &Phi1)
computes the decay angles for Vector Boson Fusion (VBF) -> H -> ZZ as defined in Figure 1 of arXiv:12...
complex(8), parameter, public c5
double Eval(const MelaIO *RcdME, TVar::VerbosityLevel verbosity) const
TVar::VerbosityLevel getVerbosity()
Gets the current verbosity level for MELA.
RooAbsReal * Lambda_z4qsq[SIZE_HVV_CQSQ]
double selfDHg4g4coupl[nSupportedHiggses][SIZE_HGG][2]
void computeP(float &prob, bool useConstant=true)
Computes the probability for the probabilities on the decay side of things using the constituent daug...
int configure(TVar::Process model_)
RooAbsReal * Lambda_z1qsq[SIZE_HVV_CQSQ]
MelaPConstant * pAvgSmooth_MCFM_Had_WH_bkgZZ_4e
std::string GetMELAPath()
void computePM4l(TVar::SuperMelaSyst syst, float &prob)
void set_InputEvent(SimpleParticleCollection_t *pDaughters, SimpleParticleCollection_t *pAssociated=0, SimpleParticleCollection_t *pMothers=0, bool isGen=false)
void setCandidateDecayMode(TVar::CandidateDecayMode mode)
Sets the decay mode for your event.
void set_CurrentCandidateFromIndex(unsigned int icand)
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_S_HSMHiggs_2mu2e
void setMelaHiggsMass(double myHiggsMass, int index=0)
Sets the mass of your chosen Higgs.
TensorPdfFactory_ppHVV * spin2Model
int getCurrentCandidateIndex()
Returns the index of the current MELA candidate - returns -1 if there is no candidate to be found.
MelaPConstant * pAvgSmooth_JHUGen_ZZGG_HSMHiggs_4mu
void makeParamsConst(bool yesNo=true)
double selfDZqqcoupl[SIZE_ZQQ][2]
RooAbsReal * ggsgs4List[1][2]
double selfDHb4b4coupl[nSupportedHiggses][SIZE_HQQ][2]
double selfDHwpwpcoupl[SIZE_HVV][2]
void computeVHAngles(float &costhetastar, float &costheta1, float &costheta2, float &Phi, float &Phi1, float &m1, float &m2, TLorentzVector p4M11, int Z1_lept1Id, TLorentzVector p4M12, int Z1_lept2Id, TLorentzVector p4M21, int Z2_lept1Id, TLorentzVector p4M22, int Z2_lept2Id, TLorentzVector jet1, int jet1Id, TLorentzVector jet2, int jet2Id, TLorentzVector *injet1=0, int injet1Id=0, TLorentzVector *injet2=0, int injet2Id=0)
double GetSigShapeSystematic(std::string parName)
MelaPConstant * pAvgSmooth_JHUGen_JJQCD_HSMHiggs[TVar::nFermionMassRemovalSchemes-1]
void computeProdXS_VH(float &mevalue, bool includeHiggsDecay=false)
RooAbsReal * gzgs2List[1][2]
MelaPConstant * pAvgSmooth_MCFM_JJVBF_S_HSMHiggs_4e
double selfDHzzLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ]
void computeProdXS_JJH(float &mevalue)
void getXPropagator(TVar::ResonancePropagatorScheme scheme, float &prop)
void computeProdDecP(double selfDHvvcoupl_input[nSupportedHiggses][SIZE_HVV][2], double selfDHwwcoupl_input[nSupportedHiggses][SIZE_HVV][2], double selfDaTQGCcoupl_input[SIZE_ATQGC][2], double selfDAZffcoupl_input[SIZE_AZff][2], float &prob, bool useConstant=true)
computes the combined production and decay probability while taking in coupling arrays
virtual void resetHypotheses()
void computeAngles(float &costhetastar, float &costheta1, float &costheta2, float &Phi, float &Phi1, TLorentzVector Z1_lept1, int Z1_lept1Id, TLorentzVector Z1_lept2, int Z1_lept2Id, TLorentzVector Z2_lept1, int Z2_lept1Id, TLorentzVector Z2_lept2, int Z2_lept2Id)
double selfDHwwcoupl[nSupportedHiggses][SIZE_HVV][2]
MelaPConstant * pAvgSmooth_MCFM_JJVBF_bkgZZ_4mu
void computeProdXS_VVHVV(float &mevalue)
void constructDggr(float bkg_VAMCFM_noscale, float ggzz_VAMCFM_noscale, float ggHZZ_prob_pure_noscale, float ggHZZ_prob_int_noscale, float widthScale, float &myDggr)
MelaPConstant * pAvgSmooth_MCFM_ZZGG_bkgZZ_4e
double GetSigShapeParameter(std::string parName)
void setRemoveJetMasses(bool MasslessLeptonSwitch=true)
either permits or forbids massive jets.
SuperDijetMela * superDijet
double selfDGvvpcoupl[SIZE_GVV][2]
void resetMass(double inmass, int ipart)
Resets the mass for a particle that is an electroweak parameter according to its id.
void set_RenFacScaleMode(TVar::EventScaleScheme renormalizationSch, TVar::EventScaleScheme factorizationSch, double ren_sf, double fac_sf)
void resetInputEvent()
Resets the event in preparation for the next iteration of the event loop.
double selfDAZffcoupl[SIZE_AZff][2]
int selfDHzzCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ]
bool configureAnalyticalPDFs()
double get_PrimaryHiggsMass()
MelaPConstant * pAvgSmooth_MCFM_JJVBF_bkgZZ_4e
void setSpinOneCouplings()
TVar::Production myProduction_
void setMelaHiggsMassWidth(double myHiggsMass, double myHiggsWidth, int index)
a combination of setMelaHiggsMass and setMelaHiggsWidth.
MELAParticle * getAssociatedJet(int index) const
void computeXS(float &mevalue)
double selfDZvvcoupl[SIZE_ZVV][2]
MelaPConstant * pAvgSmooth_JHUGen_JJVBF_HSMHiggs[TVar::nFermionMassRemovalSchemes-1]
void computeVBFAngles_ComplexBoost(float &Q2V1, float &Q2V2, float &costheta1_real, float &costheta1_imag, float &costheta2_real, float &costheta2_imag, float &Phi, float &costhetastar, float &Phi1)
virtual void makeParamsConst(bool yesNo)
MELAOutputStreamer MELAerr
ResonancePropagatorScheme
TRandom3 melaRandomNumber
void computeProdP_ttH(float &prob, int topProcess=2, int topDecay=0, bool useConstant=true)
void computeD_gg(TVar::MatrixElement myME, TVar::Process myType, float &prob)
void computeConstant(float &prob)
void writeCentered(const T &val, char fillch=' ', std::streamsize gapsize=0)
double selfDHt4t4coupl[nSupportedHiggses][SIZE_HQQ][2]
double selfDHwwLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ]
void setMelaPrimaryHiggsMass(double myHiggsMass)
Sets the mass of the "primary" higgs.
MelaPConstant * pAvgSmooth_MCFM_ZZGG_HSMHiggs_2mu2e
void computeFakeJet(TLorentzVector const &realJet, TLorentzVector const &others, TLorentzVector &fakeJet)
void setMEConst(const double &val)
void computeDijetConvBW(float &prob, bool useTrueBW=false)
complex(8), parameter, public c1
This is the "MELA" object that interfaces with the Fortran code in both MCFM-JHUGen and pure JHUGen.
RooRealVar * costhetastar_rrv
MelaPConstant * pAvgSmooth_MCFM_JJQCD_bkgZZ_4mu
std::vector< MELAParticle * > & getAssociatedJets()
RooRealVar * costheta2_rrv
double selfDHwwpcoupl[SIZE_HVV][2]
TVar::VerbosityLevel myVerbosity_
RooAbsReal * Lambda_z3qsq[SIZE_HVV_CQSQ]
MelaPConstant * pAvgSmooth_MCFM_ZZGG_HSMHiggs_4mu
RooAbsReal * gzgs1List[1][2]
bool isAGluon(const int id)
RooAbsReal * g3List[8][2]
double selfDHzzcoupl[nSupportedHiggses][SIZE_HVV][2]
std::vector< TLorentzVector > calculate4Momentum(double Mx, double M1, double M2, double theta, double theta1, double theta2, double Phi1, double Phi)
void computeP_selfDspin1(double selfDZqqcoupl_input[SIZE_ZQQ][2], double selfDZvvcoupl_input[SIZE_ZVV][2], float &prob, bool useConstant=true)
This is a function that calls Mela::computeP with preset quark, gluon, and vector boson couplings for...
void set_CurrentCandidate(MELACandidate *cand)
RooAbsReal * g1List[8][2]
TVar::FermionMassRemoval JetMassScheme
MelaPConstant * pAvgSmooth_MCFM_JJQCD_bkgZJets_2l2q
RooSpinZero::modelCouplings couplings
MelaPConstant * pAvgSmooth_MCFM_Had_ZH_bkgZZ_4mu
RooAbsReal * ggsgs3List[1][2]
VectorPdfFactory * spin1Model
double selfDHzzpcoupl[SIZE_HVV][2]
void set_LeptonInterference(TVar::LeptonInterference myLepInterf)
MelaPConstant * pAvgSmooth_MCFM_Had_WH_S_HSMHiggs_4e
void set_aTQGCCouplings(double selfDaTQGCcoupl[SIZE_ATQGC][2])
void setSpinTwoCouplings()
std::vector< SimpleParticle_t > SimpleParticleCollection_t
double selfDHggcoupl[nSupportedHiggses][SIZE_HGG][2]
RooAbsReal * Lambda_z2qsq[SIZE_HVV_CQSQ]
void reset_MCFM_EWKParameters(double ext_Gf, double ext_aemmz, double ext_mW, double ext_mZ, double ext_xW, int ext_ewscheme=3)
void setMelaLeptonInterference(TVar::LeptonInterference myLepInterf=TVar::DefaultLeptonInterf)
Sets the MELA Lepton Interference.
RooAbsReal * g4List[8][2]
void computeProdXS_ttH(float &mevalue, int topProcess, int topDecay=0)
MelaPConstant * pAvgSmooth_MCFM_JJVBF_S_HSMHiggs_2mu2e
void set_SpinTwoContact(double selfDGvvpcoupl[SIZE_GVV][2], double selfDGvpvpcoupl[SIZE_GVV][2])
double get_PrimaryWidth(int ipart)
void appendTopCandidate(SimpleParticleCollection_t *TopDaughters)
Adds a top quark as a MELA candidate.
RooAbsReal * cLambda_qsq[SIZE_HVV_CQSQ]
void SetDecayChannel(std::string myChan)
double selfDHbbcoupl[nSupportedHiggses][SIZE_HQQ][2]
TVar::LeptonInterference myLepInterf_
void computeVBFAngles(float &costhetastar, float &costheta1, float &costheta2, float &Phi, float &Phi1, float &Q2V1, float &Q2V2, TLorentzVector p4M11, int Z1_lept1Id, TLorentzVector p4M12, int Z1_lept2Id, TLorentzVector p4M21, int Z2_lept1Id, TLorentzVector p4M22, int Z2_lept2Id, TLorentzVector jet1, int jet1Id, TLorentzVector jet2, int jet2Id, TLorentzVector *injet1=0, int injet1Id=0, TLorentzVector *injet2=0, int injet2Id=0)
@ kUseAssociated_StableTops
void GetBoostedParticleVectors(MELACandidate *melaCand, TVar::simple_event_record &mela_event, TVar::VerbosityLevel verbosity=TVar::DEBUG)
RooAbsReal * bList[SIZE_GVV][2]
double selfDGggcoupl[SIZE_GGG][2]
RooAbsReal * ggsgs2List[1][2]
MelaPConstant * pAvgSmooth_MCFM_JJVBF_bkgZZ_2mu2e
MELACandidate * getCurrentCandidate()
Gets the current MELA top-level (input) candList object.
bool isAJet(const int id)
double getPrimaryMass(int ipart)
A function to get the current primary EW/QCD parameters from MELA.
void computeVBFAngles_ComplexBoost(float &costhetastar, float &costheta1_real, float &costheta1_imag, float &costheta2_real, float &costheta2_imag, float &Phi, float &Phi1, float &Q2V1, float &Q2V2, TLorentzVector p4M11, int Z1_lept1Id, TLorentzVector p4M12, int Z1_lept2Id, TLorentzVector p4M21, int Z2_lept1Id, TLorentzVector p4M22, int Z2_lept2Id, TLorentzVector jet1, int jet1Id, TLorentzVector jet2, int jet2Id, TLorentzVector *injet1=0, int injet1Id=0, TLorentzVector *injet2=0, int injet2Id=0)
complex(8), parameter, public c2
void set_PrimaryHiggsMass(double mh)
void set_wHiggs(double gah_, int index)
float GetConvBW(TVar::Production prod, MELACandidate *cand, bool useTrueBW)
MelaPConstant * pAvgSmooth_JHUGen_Had_ZH_HSMHiggs[TVar::nFermionMassRemovalSchemes-1]
void applyJetMassCorrection(bool flag=true)