53 double lheHT = 0, lheHTIncoming = 0;
54 unsigned int lheNj = 0, lheNb = 0, lheNc = 0, lheNuds = 0, lheNglu = 0;
58 const auto& hepeup = lheProd.
hepeup();
59 const auto& pup = hepeup.
PUP;
60 int lep = -1, lepBar = -1, nu = -1, nuBar = -1;
61 std::vector<float> vals_pt;
62 std::vector<float> vals_eta;
63 std::vector<float> vals_phi;
64 std::vector<float> vals_mass;
65 std::vector<float> vals_pz;
66 std::vector<int> vals_pid;
67 std::vector<int> vals_status;
68 std::vector<int> vals_spin;
69 alphaS = hepeup.AQCDUP;
70 for (
unsigned int i = 0,
n = pup.size();
i <
n; ++
i) {
74 TLorentzVector
p4(pup[
i][0], pup[
i][1], pup[
i][2], pup[
i][3]);
75 vals_pid.push_back(hepeup.IDUP[
i]);
76 vals_spin.push_back(hepeup.SPINUP[
i]);
77 vals_status.push_back(
status);
80 vals_eta.push_back(0);
81 vals_phi.push_back(0);
82 vals_mass.push_back(0);
83 vals_pz.push_back(
p4.Pz());
85 vals_pt.push_back(
p4.Pt());
86 vals_eta.push_back(
p4.Eta());
87 vals_phi.push_back(
p4.Phi());
88 vals_mass.push_back(
p4.M());
92 if ((
status == 1) && ((idabs == 21) || (idabs > 0 && idabs < 7))) {
101 else if (idabs == 21)
104 double pt = std::hypot(pup[
i][0], pup[
i][1]);
107 hepeup.MOTHUP[
i].first - 1,
109 int mothIdxTwo =
std::max(hepeup.MOTHUP[
i].second - 1, 0);
110 int mothStatus = hepeup.ISTUP[mothIdx];
111 int mothStatusTwo = hepeup.ISTUP[mothIdxTwo];
112 bool hasIncomingAsMother = mothStatus < 0 || mothStatusTwo < 0;
113 if (hasIncomingAsMother)
115 }
else if (idabs == 12 || idabs == 14 || idabs == 16) {
116 (hepeup.IDUP[
i] > 0 ? nu : nuBar) =
i;
117 }
else if (idabs == 11 || idabs == 13 || idabs == 15) {
118 (hepeup.IDUP[
i] > 0 ? lep : lepBar) =
i;
121 std::pair<int, int>
v(0, 0);
122 if (lep != -1 && lepBar != -1)
123 v = std::make_pair(lep, lepBar);
124 else if (lep != -1 && nuBar != -1)
125 v = std::make_pair(lep, nuBar);
126 else if (nu != -1 && lepBar != -1)
127 v = std::make_pair(nu, lepBar);
128 else if (nu != -1 && nuBar != -1)
129 v = std::make_pair(nu, nuBar);
130 if (
v.first != -1 &&
v.second != -1) {
131 lheVpt = std::hypot(pup[
v.first][0] + pup[
v.second][0], pup[
v.first][1] + pup[
v.second][1]);
134 out.addColumnValue<uint8_t>(
138 out.addColumnValue<uint8_t>(
140 out.addColumnValue<uint8_t>(
143 out.addColumnValue<
float>(
"HTIncoming",
145 "HT, scalar sum of parton pTs at LHE step, restricted to partons",
152 auto outPart = std::make_unique<nanoaod::FlatTable>(vals_pt.size(),
"LHEPart",
false);
154 outPart->addColumn<
float>(
156 outPart->addColumn<
float>(
158 outPart->addColumn<
float>(
160 outPart->addColumn<
float>(
163 outPart->addColumn<
int>(