|
JHUGen MELA
JHUGen v7.5.6, MELA v2.4.2
Matrix element calculations as used in JHUGen.
|
These are C++ macros that define named couplings in the Python code.
More...
|
#define | MAKE_COUPLING_ARR_SPIN_ZERO(arrayName, size, arrType) |
| Generates the array for spin 0 values in JHUGen and JHUGen-MCFM. More...
|
|
#define | MAKE_COUPLING_ARR_SPIN_ONETWO(arrayName, size) |
| Generates the array for spin 1 and spin 2 values in JHUGen. More...
|
|
#define | MAKE_COUPLING_REAL_IMAGINARY_SPIN_ZERO(arrayName, couplingName, couplingIndex, higgsIndex) |
| Generates the couplings for spin 0 values in JHUGen and JHUGen-MCFM. More...
|
|
#define | MAKE_COUPLING_REAL_IMAGINARY_SPIN_ONETWO(arrayName, couplingName, couplingIndex) |
| Generates the couplings for spin 1 and spin 2 values in JHUGen. More...
|
|
#define | MAKE_COUPLING_C_LAMBDA(arrayName, couplingName, couplingIndex, higgsIndex) |
| Generates the couplings for C Lambda values in JHUGen and JHUGen-MCFM. More...
|
|
#define | MAKE_COUPLING_LAMBDA(arrayName, couplingName, couplingIndex_1, couplingIndex_2, higgsIndex) |
| Generates the couplings for Lambda values in JHUGen and JHUGen-MCFM. More...
|
|
#define | MAKE_COUPLING_MADMELA(couplingName, couplingIndex_1) |
| Generates the couplings for SMEFTSim Wilson Coefficients in MadMELA. More...
|
|
These are C++ macros that define named couplings in the Python code.
◆ MAKE_COUPLING_ARR_SPIN_ONETWO
#define MAKE_COUPLING_ARR_SPIN_ONETWO |
( |
|
arrayName, |
|
|
|
size |
|
) |
| |
Value: .def(#arrayName, [](py::object &
obj){ \
return py::array_t<double>(std::vector<double>{size, 2}, (
const double*) &D.arrayName,
obj); \
})
Generates the array for spin 1 and spin 2 values in JHUGen.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
size | This is the size of the array |
- Note
- These coupling entries have both a real and an imaginary component. You can set them via:
import Mela
m.<arrayName>()[<couplingIndex>][0] = <real>
m.<arrayName>()[<couplingIndex>][1] = <imag>
Definition at line 418 of file mela_binding.cpp.
◆ MAKE_COUPLING_ARR_SPIN_ZERO
#define MAKE_COUPLING_ARR_SPIN_ZERO |
( |
|
arrayName, |
|
|
|
size, |
|
|
|
arrType |
|
) |
| |
Value: .def(#arrayName, [](py::object &
obj){ \
return py::array_t<arrType>(std::vector<arrType>{
nSupportedHiggses, size, 2}, (
const arrType*) &D.arrayName,
obj); \
})
Generates the array for spin 0 values in JHUGen and JHUGen-MCFM.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
size | This is the size of the array |
arrType | The data type of the array |
- Note
- These coupling entries have both a real and an imaginary component. You can set them via:
import Mela
m.<arrayName>()[<couplingIndex>][0] = <real>
m.<arrayName>()[<couplingIndex>][1] = <imag>
Definition at line 399 of file mela_binding.cpp.
◆ MAKE_COUPLING_C_LAMBDA
#define MAKE_COUPLING_C_LAMBDA |
( |
|
arrayName, |
|
|
|
couplingName, |
|
|
|
couplingIndex, |
|
|
|
higgsIndex |
|
) |
| |
Value: .def_property(\
#couplingName, \
py::cpp_function(\
return array_val.at(higgsIndex,couplingIndex);\
}),\
py::cpp_function(\
[](py::object &
obj,
double coupl){\
array_val.mutable_at(higgsIndex,couplingIndex) = coupl;\
}, py::keep_alive<0, 1>())\
)
Generates the couplings for C Lambda values in JHUGen and JHUGen-MCFM.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
couplingName | This is the name of the coupling to set (should be unique for every couplingIndex/higgsIndex combination) |
couplingIndex | This is the index at which the name corresponds to the array |
higgsIndex | This is the index of which Higgs to use (MCFM supports 2 resonances) |
- Note
- These coupling entries only have a real component. You can set them via:
import Mela
m.<couplingName> = value
print(m.<couplingName>)
Definition at line 498 of file mela_binding.cpp.
◆ MAKE_COUPLING_LAMBDA
#define MAKE_COUPLING_LAMBDA |
( |
|
arrayName, |
|
|
|
couplingName, |
|
|
|
couplingIndex_1, |
|
|
|
couplingIndex_2, |
|
|
|
higgsIndex |
|
) |
| |
Value: .def_property(\
#couplingName, \
py::cpp_function(\
return array_val.at(higgsIndex,couplingIndex_1, couplingIndex_2);\
}),\
py::cpp_function(\
[](py::object &
obj,
double coupl){\
array_val.mutable_at(higgsIndex,couplingIndex_1, couplingIndex_2) = coupl;\
}, py::keep_alive<0, 1>())\
)
Generates the couplings for Lambda values in JHUGen and JHUGen-MCFM.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
couplingName | This is the name of the coupling to set (should be unique for every couplingIndex/higgsIndex combination) |
couplingIndex_1 | This is the index at which the name corresponds to the array |
couplingIndex_2 | This is the index at which the name corresponds to the array |
higgsIndex | This is the index of which Higgs to use (MCFM supports 2 resonances) |
- Note
- These coupling entries only have a real component. You can set them via:
import Mela
m.<couplingName> = value
print(m.<couplingName>)
Definition at line 532 of file mela_binding.cpp.
◆ MAKE_COUPLING_MADMELA
#define MAKE_COUPLING_MADMELA |
( |
|
couplingName, |
|
|
|
couplingIndex_1 |
|
) |
| |
Value: .def_property(\
#couplingName, \
py::cpp_function(\
py::array_t array_val = py::array_t<double>(std::vector<double>{
SIZE_SMEFT}, (
const double*) &D.selfDSMEFTSimcoupl,
obj);\
return array_val.at(couplingIndex_1);\
}),\
py::cpp_function(\
[](py::object &
obj,
double coupl){\
py::array_t array_val = py::array_t<double>(std::vector<double>{
SIZE_SMEFT}, (
const double*) &D.selfDSMEFTSimcoupl,
obj);\
array_val.mutable_at(couplingIndex_1) = coupl;\
}, py::keep_alive<0, 1>())\
)
Generates the couplings for SMEFTSim Wilson Coefficients in MadMELA.
- Parameters
-
couplingName | This is the name of the Wilson Coefficient to set |
couplingIndex_1 | This is the index at which the name corresponds to the array |
- Note
- These coupling entries only have a real component. You can set them via:
import Mela
m.<couplingName> = value
print(m.<couplingName>)
Definition at line 562 of file mela_binding.cpp.
◆ MAKE_COUPLING_REAL_IMAGINARY_SPIN_ONETWO
#define MAKE_COUPLING_REAL_IMAGINARY_SPIN_ONETWO |
( |
|
arrayName, |
|
|
|
couplingName, |
|
|
|
couplingIndex |
|
) |
| |
Value: .def_property(\
#couplingName,\
py::cpp_function(\
return py::array_t<double>(std::vector<double>{2}, (
const double*) &D.arrayName[couplingIndex],
obj);\
}, py::keep_alive<0, 1>()),\
py::cpp_function(\
[](
Mela &D, std::array<double, 2> coupl){\
D.arrayName[couplingIndex][0] = coupl[0];\
D.arrayName[couplingIndex][1] = coupl[1];\
}, py::keep_alive<0, 1>())\
)
Generates the couplings for spin 1 and spin 2 values in JHUGen.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
couplingName | This is the name of the coupling to set (should be unique for every couplingIndex/higgsIndex combination) |
couplingIndex | This is the index at which the name corresponds to the array |
- Note
- These coupling entries have both a real and an imaginary component. You can set them via:
import Mela
m.<couplingName> = [real, imag]
print(m.<couplingName>)
Definition at line 468 of file mela_binding.cpp.
◆ MAKE_COUPLING_REAL_IMAGINARY_SPIN_ZERO
#define MAKE_COUPLING_REAL_IMAGINARY_SPIN_ZERO |
( |
|
arrayName, |
|
|
|
couplingName, |
|
|
|
couplingIndex, |
|
|
|
higgsIndex |
|
) |
| |
Value: .def_property(\
#couplingName,\
py::cpp_function(\
return py::array_t<double>(std::vector<double>{2}, (
const double*) &D.arrayName[higgsIndex][couplingIndex],
obj);\
}, py::keep_alive<0, 1>()),\
py::cpp_function(\
[](
Mela &D, std::array<double, 2> coupl){\
D.arrayName[higgsIndex][couplingIndex][0] = coupl[0];\
D.arrayName[higgsIndex][couplingIndex][1] = coupl[1];\
}, py::keep_alive<0, 1>())\
)
Generates the couplings for spin 0 values in JHUGen and JHUGen-MCFM.
- Parameters
-
arrayName | This is the name of the array containing the coupling |
couplingName | This is the name of the coupling to set (should be unique for every couplingIndex/higgsIndex combination) |
couplingIndex | This is the index at which the name corresponds to the array |
higgsIndex | This is the index of which Higgs to use (MCFM supports 2 resonances) |
- Note
- These coupling entries have both a real and an imaginary component. You can set them via:
import Mela
m.<couplingName> = [real, imag]
print(m.<couplingName>)
Definition at line 439 of file mela_binding.cpp.
◆ PYBIND11_MODULE()
PYBIND11_MODULE |
( |
Mela |
, |
|
|
m |
|
|
) |
| |
The actual binding code for MELA.
Definition at line 583 of file mela_binding.cpp.
584 py::class_<SimpleParticle_t>(m,
"SimpleParticle_t")
585 .def(py::init(&
particle_initializer), py::arg(
"id"), py::arg(
"x"), py::arg(
"y"), py::arg(
"z"), py::arg(
"e"), py::arg(
"ptEtaPhi") =
false)
594 return py::make_tuple(P.second.Px(), P.second.Py(), P.second.Pz(), P.second.E());
597 return py::make_tuple(P.second.Pt(), P.second.Eta(), P.second.Phi(), P.second.M());
600 P.second = TLorentzVector(Px, Py, Pz, E);
603 return "SimpleParticle(id=" + std::to_string(P.first) +
",P4=<" + std::to_string(P.second.Px()) +
", " + std::to_string(P.second.Py()) +
", " + std::to_string(P.second.Pz()) +
", " + std::to_string(P.second.E()) +
">)";
607 py::class_<SimpleParticleCollection_t>(m,
"SimpleParticleCollection_t")
615 return py::make_iterator(C.begin(), C.end());
616 }, py::keep_alive<0, 1>())
618 py::list list_type = py::cast(C);
628 TLorentzVector sum = TLorentzVector(0,0,0,0);
632 return py::make_tuple(sum.Px(), sum.Py(), sum.Pz(), sum.E());
635 TLorentzVector sum = TLorentzVector(0,0,0,0);
644 for(
int i = 0; i < (int)C.size(); i++){
645 pickleable.append(py::make_tuple(C[i].first, C[i].
second.Px(), C[i].second.Py(), C[i].second.Pz(), C[i].second.E()));
647 return py::cast<py::tuple>(pickleable);
650 std::vector<int> ids;
651 std::vector<double> x;
652 std::vector<double> y;
653 std::vector<double>
z;
654 std::vector<double> e;
656 for(
int i = 0; i < (int)t.size(); i++){
658 ids.push_back(o[0].cast<int>());
659 x.push_back(o[1].cast<double>());
660 y.push_back(o[2].cast<double>());
661 z.push_back(o[3].cast<double>());
662 e.push_back(o[4].cast<double>());
668 return "SimpleParticleCollection_t(length=" + to_string(C.size()) +
")";
671 py::class_<MELAParticle>(m,
"MELAParticle")
673 .def(py::init<int>())
674 .def(py::init([](
int id_,
double Px,
double Py,
double Pz,
double E){
675 TLorentzVector vec = TLorentzVector();
676 vec.SetPxPyPzE(Px, Py, Pz, E);
679 .def(py::init<MELAParticle const &>())
680 .
def(
"assign", &MELAParticle::operator=)
684 .
def(
"addVec", [](
MELAParticle& mp,
double Px,
double Py,
double Pz,
double E){
685 TLorentzVector vec = TLorentzVector();
686 vec.SetPxPyPzE(Px, Py, Pz, E);
690 .def_property_readonly(
"p4",
692 TLorentzVector p4 = mp.
p4;
693 return py::make_tuple(p4.Px(), p4.Py(), p4.Pz(), p4.E());
695 .
def(
"setP4", [](
MELAParticle& mp,
double Px,
double Py,
double Pz,
double E){
696 TLorentzVector vec = TLorentzVector();
697 vec.SetPxPyPzE(Px, Py, Pz, E);
729 .def(
"dot", [](
MELAParticle& m,
double Px,
double Py,
double Pz,
double E){
730 TLorentzVector vec = TLorentzVector();
731 vec.SetPxPyPzE(Px, Py, Pz, E);
737 .
def(
"euclidean_dot", [](
MELAParticle& m,
double Px,
double Py,
double Pz,
double E){
738 TLorentzVector vec = TLorentzVector();
739 vec.SetPxPyPzE(Px, Py, Pz, E);
745 .
def(
"deltaR", [](
MELAParticle& m,
double Px,
double Py,
double Pz,
double E){
746 TLorentzVector vec = TLorentzVector();
747 vec.SetPxPyPzE(Px, Py, Pz, E);
753 .
def(
"boost", [](
MELAParticle& m,
double Px,
double Py,
double Pz,
bool boostAll){
754 TVector3 v = TVector3(Px, Py, Pz);
755 m.
boost(v, boostAll);
759 TVector3 v = m.
vect();
760 return py::make_tuple(v.X(), v.Y(), v.Z());
764 return py::make_tuple(v.X(), v.Y(), v.Z());
768 return "MELAParticle(id=" + to_string(mp.
id) +
",passSelection=" + to_string(mp.
passSelection) +
",genStatus=" + to_string(mp.
genStatus) +
",lifetime=" + to_string(mp.
lifetime) +
")";
771 py::class_<MELAThreeBodyDecayCandidate, MELAParticle>(m,
"MELAThreeBodyDecayCandidate")
773 .def(py::init([](
int id_,
double Px,
double Py,
double Pz,
double E){
774 TLorentzVector vec = TLorentzVector();
775 vec.SetPxPyPzE(Px, Py, Pz, E);
778 .def(py::init<MELAParticle*, MELAParticle*, MELAParticle*>())
779 .
def(py::init<MELAThreeBodyDecayCandidate const &>())
780 .
def(
"assign", &MELAThreeBodyDecayCandidate::operator=)
786 return mt.getPartnerParticle();
789 return mt.getWFermion();
792 return mt.getWAntifermion();
798 py::class_<MELACandidate, MELAParticle>(m,
"MELACandidate")
800 .def(py::init<int, bool>(), py::arg(
"id_"), py::arg(
"associatedByHighestPt_")=
false)
801 .def(py::init([](
int id_,
double Px,
double Py,
double Pz,
double E,
bool associatedByHighestPt_){
802 TLorentzVector vec = TLorentzVector();
803 vec.SetPxPyPzE(Px, Py, Pz, E);
805 }), py::arg(
"id_"), py::arg(
"Px"), py::arg(
"Py"), py::arg(
"Pz"), py::arg(
"E"), py::arg(
"associatedByHighestPt_")=
false)
806 .
def(py::init<MELACandidate const &>())
807 .
def(
"assign", &MELACandidate::operator=)
818 return mc.getSortedDaughters();
821 return mc.getSortedVs();
824 return mc.getAssociatedLeptons();
827 return mc.getAssociatedNeutrinos();
830 return mc.getAssociatedPhotons();
833 return mc.getAssociatedJets();
836 return mc.getAssociatedTops();
839 return mc.getAssociatedSortedVs();
866 mc.addUnordered(myParticle, particleArray);
869 mc.addUnordered(myParticle, particleArray);
872 mc.addByHighestPt(myParticle, particleArray);
875 mc.addByHighestPt(myParticle, particleArray);
878 py::class_<TVar::event_scales_type>(m,
"event_scales_type")
879 .def(py::init<TVar::EventScaleScheme, TVar::EventScaleScheme, double, double>())
885 py::class_<TVar::simple_event_record>(m,
"simple_event_record")
912 py::class_<Mela>(m,
"Mela")
913 .def(py::init<double, double, TVar::VerbosityLevel>())
914 .def(py::init<double, double>())
915 .def(py::init<double>())
917 .def(py::init<Mela const &>())
918 .def(
"__repr__", [](
Mela&
me){
923 .def(
"setInputEvent", &
Mela::setInputEvent, py::arg(
"pDaughters"), py::arg(
"pAssociated")=
nullptr, py::arg(
"pMothers")=
nullptr, py::arg(
"isGen")=
false)
930 .def(
"SetMadgraphCKMElements", &
Mela::SetMadgraphCKMElements, py::arg(
"ckmlambda")=0.2265, py::arg(
"ckma")=0.79, py::arg(
"ckmrho")=0.141, py::arg(
"ckmeta")=0.357, py::arg(
"force_refresh")=
false)
952 .def(
"computeP", &
computeP, py::arg(
"useConstant"))
965 .def(
"PrintCurrentCandidateSummary", [](
Mela& mela){
989 .def(
"selfDHzzLambda_qsq", [](py::object &
obj){ \
993 .
def(
"selfDHwwLambda_qsq", [](py::object &
obj){ \
997 .
def(
"selfDHzzCLambda_qsq", [](py::object &
obj){ \
1001 .
def(
"selfDHwwCLambda_qsq", [](py::object &
obj){ \
1002 Mela &D =
obj.cast<
Mela&>(); \
1020 .def(
"selfDSMEFTSimcoupl", [](py::object &
obj){ \
1021 Mela &D =
obj.cast<
Mela&>(); \
1022 return py::array_t<double>(std::vector<double>{
SIZE_SMEFT}, (
const double*) &D.selfDSMEFTSimcoupl,
obj); \
complex(8), public ghzpzp1_prime7
float getXPropagator(Mela &mela, TVar::ResonancePropagatorScheme scheme)
analog of Mela::getXPropagator
complex(8), public bzpgs8
complex(8), public ghzpzp4_prime4
std::pair< int, TLorentzVector > SimpleParticle_t
double selfDHg4g4coupl[nSupportedHiggses][SIZE_HGG][2]
MELATopCandidate_t * getAssociatedTop(int index) const
complex(8), public ghzzp1_prime
complex(8), public bzpgs3
py::array getWeightedMEArray(Mela &mela)
the analog of MelaIO::getWeightedMEArray
bool daughtersInterfere() const
std::vector< SimpleParticleCollection_t > pTopDaughters
SimpleParticleCollection_t pStableAntitops
complex(8), public ghzpzp2_prime7
pymela::gHIGGS_KAPPA value("gHIGGS_KAPPA_TILDE", pymela::gHIGGS_KAPPA_TILDE) .value("SIZE_HQQ"
MELAParticle * getDaughter(int index) const
complex(8), public ghzzp1_prime7
void setAddAssociatedByHighestPt(bool associatedByHighestPt_)
complex(8), public ghzpzp1
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.
double selfDHzzpcoupl[SIZE_HVV][2]
std::complex< double > GetMadgraphCKMElement(int iquark, int jquark)
float compute4FermionWeight(Mela &mela)
analog of Mela::compute4FermionWeight
complex(8), public bzpzp1
double selfDaTQGCcoupl[SIZE_ATQGC][2]
complex(8), public ghzzp3_prime
complex(8), public bzpzp6
complex(8), public ghzzp4_prime2
double selfDHzzLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ]
double selfDGggcoupl[SIZE_GGG][2]
void addAssociatedJet(MELAParticle *myParticle)
float computeD_CP(Mela &mela, TVar::MatrixElement myME, TVar::Process myType)
analog of Mela::computeD_CP
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.
complex(8), public ghwpwp1
double euclidean_dot(const TLorentzVector &v) const
double selfDZpffcoupl[SIZE_Vpff][2]
int getNDaughters() const
double selfDHb4b4coupl[nSupportedHiggses][SIZE_HQQ][2]
complex(8), public zprime_zz_1
int nRequested_AssociatedJets
double selfDZqqcoupl[SIZE_ZQQ][2]
complex(8), public bzpgs2
void swap(MELAParticle &particle_)
complex(8), public ghzpzp2_prime6
complex(8), public ghzzp2_prime5
@ nFermionMassRemovalSchemes
complex(8), public ghzpzp2
complex(8), public ghzzp4_prime3
complex(8), public bzpzp8
#define MAKE_COUPLING_REAL_IMAGINARY_SPIN_ZERO(arrayName, couplingName, couplingIndex, higgsIndex)
Generates the couplings for spin 0 values in JHUGen and JHUGen-MCFM.
MELAParticle * getAssociatedNeutrino(int index) const
complex(8), public ghzzp4_prime4
MelaIO * getIORecord()
Returns the MELAIO object, and by consequence, the entire parton-by-parton matrix element record.
complex(8), public ghzpzp3_prime5
complex(8), public ghzpzp3_prime2
MELAParticle * getSortedV(int index) const
static bool checkCandidateExists(MELAThreeBodyDecayCandidate const *myParticle, std::vector< MELAThreeBodyDecayCandidate * > const &particleArray)
double selfDGqqcoupl[SIZE_GQQ][2]
int getNAssociatedJets() const
void addAssociatedNeutrino(MELAParticle *myParticle)
void setGenStatus(int status_)
int getNAssociatedPhotons() const
void setVerbosity(TVar::VerbosityLevel verbosity_=TVar::ERROR)
Sets the verbosity for MELA outside of the initial constructor.
complex(8), public ghzpzp1_prime5
float computePM4l(Mela &mela, TVar::SuperMelaSyst syst)
analog of Mela::computePM4l
void swap(MELACandidate &particle_)
void resetWidth(double inwidth, int ipart)
Resets the width for a particle that is an electroweak parameter according to its id.
double selfDHt4t4coupl[nSupportedHiggses][SIZE_HQQ][2]
void setMelaHiggsWidth(double myHiggsWidth=-1, int index=0)
Sets the width of your chosen Higgs.
complex(8), public ghzpzp3_prime
int getNAssociatedTops() const
#define MAKE_COUPLING_LAMBDA(arrayName, couplingName, couplingIndex_1, couplingIndex_2, higgsIndex)
Generates the couplings for Lambda values in JHUGen and JHUGen-MCFM.
complex(8), public bzpzp3
void setRenFacScaleMode(TVar::EventScaleScheme renormalizationSch, TVar::EventScaleScheme factorizationSch, double ren_sf, double fac_sf)
Sets the renormalization and the factorization schemes.
SimpleParticleCollection_t pDaughters
complex(8), public bgsgs1
complex(8), public bgsgs8
complex(8), public graviton_qq_left
complex(8), public ghzzp3_prime3
complex(8), public ghzpzp1_prime4
void setWAntifermion(MELAParticle *myParticle)
#define MAKE_COUPLING_REAL_IMAGINARY_SPIN_ONETWO(arrayName, couplingName, couplingIndex)
Generates the couplings for spin 1 and spin 2 values in JHUGen.
complex(8), public ghzzp1_prime3
complex(8), public ghzpzp1_prime6
complex(8), public bzzp10
complex(8), public bgsgs4
TVector3 calculateTotalDisplacement() const
py::array getUnweightedMEArray(Mela &mela)
the analog of MelaIO::getUnweightedMEArray
double selfDHbbcoupl[nSupportedHiggses][SIZE_HQQ][2]
def("selfDHggcoupl", [](py::object &obj){ Mela &D=obj.cast< Mela & >();return py::array_t< double >(std::vector< int >{nSupportedHiggses, SIZE_HGG, 2},(const double *) &D.selfDHggcoupl, obj);}) .def("selfDHg4g4coupl"
MELAParticle * getAssociatedPhoton(int index) const
complex(8), public ghzpzp4_prime
complex(8), public ghzpzp4_prime5
static void cleanLinkedFiles()
void PrintCandidateSummary(MELACandidate *cand)
complex(8), public ghzzp1
void addAssociatedPhoton(MELAParticle *myParticle)
complex(8), public ghzzp1_prime6
complex(8), public ghzzp4
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....
int nRequested_AssociatedPhotons
double selfDHzpzpcoupl[SIZE_HVV][2]
real(8), parameter, public second
Mela(double LHCsqrts_=13., double mh_=125., TVar::VerbosityLevel verbosity_=TVar::ERROR)
the MELA constructor
void resetQuarkMasses()
Resets the masses of each quark to their original values.
complex(8), public zprime_qq_right
complex(8), public ghzpgs3
int AssociationVCompatibility
int selfDHwwCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ]
void getRelatedParticles(std::vector< MELAParticle * > &particles) const
py::tuple getPartonWeights(Mela &mela)
the analog of MelaIO::getPartonWeights
void swap(MELAThreeBodyDecayCandidate &particle_)
complex(8), public bzpzp9
MELACandidate * shallowCopy()
SimpleParticleCollection_t collection_initializer(py::list listOfParticles)
This function initializes a single SimpleParticleCollection_t (or a list of SimpleParticle_t) in the ...
complex(8), public ghzzp2_prime6
SimpleParticle_t particle_initializer(int id, float x, float y, float z, float e, bool ptEtaPhi=false)
This function intializes a single SimpleParticle_t in the Python.
void setWFermion(MELAParticle *myParticle)
SimpleParticleCollection_t pMothers
double selfDHzzcoupl[nSupportedHiggses][SIZE_HVV][2]
array< float, 8 > computeDecayAngles(Mela &mela)
the analog of Mela::computeDecayAngles
complex(8), public ghzzp4_prime5
complex(8), public bzpzp5
complex(8), public ghzzp2_prime
bool differentiate_HWW_HZZ
float computeD_gg(Mela &mela, TVar::MatrixElement myME, TVar::Process myType)
analog of Mela::computeD_gg
complex(8), public ghzzp2
double getHiggsWidthAtPoleMass(double mass)
Returns the width of the Higgs at a given pole mass as a calculation.
complex(8), public ghzzp3_prime7
complex(8), public ghwwp1
complex(8), public ghzzp1_prime4
std::vector< int > intermediateVid
complex(8), public ghzpgs4
complex(8), public ghzpzp2_prime3
complex(8), public zprime_qq_left
double selfDHqqcoupl[nSupportedHiggses][SIZE_HQQ][2]
complex(8), public bzpgs4
float computeDijetConvBW(Mela &mela, bool useTrueBW)
the analog of Mela::computeDijetConvBW
complex(8), public ghzpzp3
void addSortedV(MELAParticle *myParticle)
SimpleParticleCollection_t pStableTops
complex(8), public ghzpzp4_prime3
float getPAux(Mela &mela)
the analog of Mela::getPAux
MELAParticle * getSortedDaughter(int index) const
MELAParticle * getAssociatedLepton(int index) const
complex(8), public ghzpzp4_prime6
void setCandidateDecayMode(TVar::CandidateDecayMode mode)
Sets the decay mode for your event.
double deltaR(const TLorentzVector &v) const
SimpleParticleCollection_t collection_initializer_from_column(std::vector< int > ids, std::vector< double > x, std::vector< double > y, std::vector< double > z, std::vector< double > e, bool ptEtaPhi=false)
This function initializes a single SimpleParticleCollection_t (or a list of SimpleParticle_t) in the ...
void setMelaHiggsMass(double myHiggsMass, int index=0)
Sets the mass of your chosen Higgs.
complex(8), public ghzzp2_prime3
void setLifetime(double life_)
complex(8), public ghzpzp2_prime2
complex(8), public ghzpzp1_prime
complex(8), public ghzzp3
int selfDHzzCLambda_qsq[nSupportedHiggses][SIZE_HVV_CQSQ]
void setPartnerParticle(MELAParticle *myParticle)
CouplingIndex_HVV
This is the enumerator for the couplings between the Higgs and the vector bosons (Z/Z' & W/W')
complex(8), public bzpzp4
void setDecayMode(TVar::CandidateDecayMode flag)
complex(8), public ghzpzp4_prime2
complex(8), public kappa_tilde
array< float, 7 > computeVHAngles(Mela &mela, TVar::Production prod)
analog of Mela::computeVHAngles
CouplingIndex_HQQ
This is the enumeration for couplings between the Higgs and the kappa formulation of quarks.
complex(8), public ghzzp3_prime5
complex(8), public ghzpgs2
CouplingIndex_HGG
This is the enumerator for the couplings between the Higgs and gluons.
void SetMadgraphCKMElements(double ckmlambda=0.2265, double ckma=0.79, double ckmrho=0.141, double ckmeta=0.357, bool force_refresh=false)
double selfDHggcoupl[nSupportedHiggses][SIZE_HGG][2]
float computeProdP_ttH(Mela &mela, int topProcess, int topDecay, bool useConstant)
analog of Mela::computeProdP_ttH
virtual void getDaughterParticles(std::vector< MELAParticle * > &particles) const
complex(8), public bzpzp7
void getDaughterParticles(std::vector< MELAParticle * > &particles) const
float computeProdP(Mela &mela, bool useConstant=true)
analog of Mela::computeProdP
complex(8), public ghzzp2_prime7
complex(8), public ghzpzp3_prime7
void resetMass(double inmass, int ipart)
Resets the mass for a particle that is an electroweak parameter according to its id.
void addAssociatedTop(MELATopCandidate_t *myParticle)
complex(8), public ghzpzp2_prime5
void resetInputEvent()
Resets the event in preparation for the next iteration of the event loop.
TVar::EventScaleScheme factorizationScheme
complex(8), public ghzzp4_prime7
void addAssociatedLepton(MELAParticle *myParticle)
int getNAssociatedNeutrinos() const
void setMelaHiggsMassWidth(double myHiggsMass, double myHiggsWidth, int index)
a combination of setMelaHiggsMass and setMelaHiggsWidth.
double selfDHwwpcoupl[SIZE_HVV][2]
MELAParticle * getAssociatedJet(int index) const
complex(8), public ghzpzp4_prime7
void addMother(MELAParticle *myParticle)
complex(8), public ghzzp1_prime2
complex(8), public bzpzp10
void resetYukawaMass(double inmass, int ipart)
Resets the Yukawa coupling (or "Yukawa mass") for a specific particle. Does not change its "intrinsic...
This is the main Mela object!
CouplingIndex_LAMBDAQSQ
This is the enumeration for couplings between the Higgs and the lambda Z/W couplings.
#define MAKE_COUPLING_ARR_SPIN_ZERO(arrayName, size, arrType)
Generates the array for spin 0 values in JHUGen and JHUGen-MCFM.
float computeProdP_VH(Mela &mela, bool includeHiggsDecay, bool useConstant)
analog of Mela::computeProdP_VH
complex(8), public graviton_qq_right
double selfDHttcoupl[nSupportedHiggses][SIZE_HQQ][2]
ResonancePropagatorScheme
const TVar::event_scales_type & getRenFacScaleMode() const
Returns the current RenFac Scale Mode for MELA.
complex(8), public bzpzp2
array< float, 7 > computeVBFAngles(Mela &mela)
analog of Mela::computeVBFAngles
complex(8), public ghzpzp1_prime2
complex(8), public ghzzp3_prime2
complex(8), public ghzpzp2_prime4
complex(8), public ghzpzp1_prime3
complex(8), public ghzpgs1_prime2
complex(8), public ghzpzp3_prime4
#define MAKE_COUPLING_C_LAMBDA(arrayName, couplingName, couplingIndex, higgsIndex)
Generates the couplings for C Lambda values in JHUGen and JHUGen-MCFM.
complex(8), public ghzzp2_prime4
SimpleParticleCollection_t pAssociated
complex(8), public ghzzp3_prime6
complex(8), public ghzpzp2_prime
double selfDWpffcoupl[SIZE_Vpff][2]
std::vector< TLorentzVector > calculate4Momentum(double Mx, double M1, double M2, double theta, double theta1, double theta2, double Phi1, double Phi)
MELAParticle * getMother(int index) const
float getConstant(Mela &mela)
the analog of Mela::getConstant
complex(8), public bzpgs1
double selfDGvpvpcoupl[SIZE_GVV][2]
TVar::EventScaleScheme renomalizationScheme
double selfDAZffcoupl[SIZE_AZff][2]
complex(8), public bgsgs2
std::vector< SimpleParticleCollection_t > pAntitopDaughters
float computeP(Mela &mela, bool useConstant=true)
analog of Mela::computeP
complex(8), public ghzpzp3_prime6
float computeProdDecP(Mela &mela, bool useConstant=true)
analog of MelThey cover a slew of uses, from converting a::computeProdDecP
int nRequested_AssociatedLeptons
double complex, dimension(2, 2) z
std::vector< SimpleParticle_t > SimpleParticleCollection_t
void testPreSelectedDaughters()
double dot(const TLorentzVector &v) const
double selfDHwpwpcoupl[SIZE_HVV][2]
virtual void getRelatedParticles(std::vector< MELAParticle * > &particles) const
complex(8), public zprime_zz_2
complex(8), public ghzzp1_prime5
void setMelaLeptonInterference(TVar::LeptonInterference myLepInterf=TVar::DefaultLeptonInterf)
Sets the MELA Lepton Interference.
complex(8), public ghzpzp3_prime3
complex(8), public ghzpzp4
double selfDHwwcoupl[nSupportedHiggses][SIZE_HVV][2]
complex(8), public ghzzp2_prime2
#define MAKE_COUPLING_ARR_SPIN_ONETWO(arrayName, size)
Generates the array for spin 1 and spin 2 values in JHUGen.
complex(8), public ghzzp3_prime4
complex(8), public ghzzp4_prime6
double selfDGvvpcoupl[SIZE_GVV][2]
array< float, 9 > computeVBFAngles_ComplexBoost(Mela &mela)
analog of Mela::computeVBFAngles_ComplexBoost
int getNAssociatedLeptons() const
virtual std::vector< int > getDaughterIds() const
double selfDHwwLambda_qsq[nSupportedHiggses][SIZE_HVV_LAMBDAQSQ][SIZE_HVV_CQSQ]
void boost(const TVector3 &vec, bool boostAll=false)
double selfDZvvcoupl[SIZE_ZVV][2]
complex(8), public bgsgs3
void setSelected(bool isSelected=true)
@ Dynamic_RandomUniform_Constrained
MELACandidate * getCurrentCandidate()
Gets the current MELA top-level (input) candList object.
#define MAKE_COUPLING_MADMELA(couplingName, couplingIndex_1)
Generates the couplings for SMEFTSim Wilson Coefficients in MadMELA.
double getPrimaryMass(int ipart)
A function to get the current primary EW/QCD parameters from MELA.
double selfDGvvcoupl[SIZE_GVV][2]
void testPreSelectedDaughters()
complex(8), public ghzzp4_prime
void setShallowCopy(bool flag)
void addDaughter(MELAParticle *myParticle)