20 #include "fastjet/PseudoJet.hh"
21 #include "fastjet/ClusterSequence.hh"
46 ATH_CHECK(m_dressParticlesKey.initialize());
47 m_decorator_eKey = m_dressParticlesKey.key() +
".e_dressed";
49 m_decorator_ptKey = m_dressParticlesKey.key() +
".pt_dressed";
50 ATH_CHECK(m_decorator_ptKey.initialize());
51 m_decorator_etaKey= m_dressParticlesKey.key() +
".eta_dressed";
52 ATH_CHECK(m_decorator_etaKey.initialize());
53 m_decorator_phiKey= m_dressParticlesKey.key() +
".phi_dressed";
54 ATH_CHECK(m_decorator_phiKey.initialize());
55 m_decorator_pt_visKey = m_dressParticlesKey.key() +
".pt_vis_dressed";
56 ATH_CHECK(m_decorator_pt_visKey.initialize());
57 m_decorator_eta_visKey= m_dressParticlesKey.key() +
".eta_vis_dressed";
58 ATH_CHECK(m_decorator_eta_visKey.initialize());
59 m_decorator_phi_visKey= m_dressParticlesKey.key() +
".phi_vis_dressed";
60 ATH_CHECK(m_decorator_phi_visKey.initialize());
61 m_decorator_m_visKey= m_dressParticlesKey.key() +
".m_vis_dressed";
62 ATH_CHECK(m_decorator_m_visKey.initialize());
63 m_decorator_nphotonKey = m_dressParticlesKey.key() +
".nPhotons_dressed";
64 ATH_CHECK(m_decorator_nphotonKey.initialize());
65 if (!m_decorationName.empty()) {m_decorationKey = m_particlesKey.key()+
"."+m_decorationName;}
66 else {m_decorationKey = m_particlesKey.key()+
".unusedPhotonDecoration";}
72 if (m_listOfPIDs.size() > 1) {
73 if (
std::find(m_listOfPIDs.begin(), m_listOfPIDs.end(), 15) != m_listOfPIDs.end()) {
74 ATH_MSG_ERROR(
"Truth taus must be dressed separately from other truth particles");
75 return StatusCode::FAILURE;
79 return StatusCode::SUCCESS;
86 const EventContext& ctx = Gaudi::Hive::currentContext();
90 if (!truthParticles.
isValid()) {
91 ATH_MSG_ERROR(
"Couldn't retrieve TruthParticle collection with name " << m_particlesKey);
92 return StatusCode::FAILURE;
96 if (!dressTruthParticles.
isValid()) {
97 ATH_MSG_ERROR(
"Couldn't retrieve TruthParticle collection with name " << m_dressParticlesKey);
98 return StatusCode::FAILURE;
115 if (!m_decorationName.empty()){
116 for (
const auto *
particle : *truthParticles){
130 std::vector<const xAOD::TruthParticle*> listOfParticlesToDress;
131 std::vector<xAOD::TruthParticle::FourMom_t> listOfDressedP4;
132 std::vector<xAOD::TruthParticle::FourMom_t> listOfBareP4;
134 if (m_listOfPIDs.size()==1 && std::abs(m_listOfPIDs[0])==15) {
136 for (
auto *pItr : *dressTruthParticles) {
137 if (!pItr->isTau()) {
139 return StatusCode::FAILURE;
142 ATH_MSG_ERROR(
"Visible momentum not available for truth taus, cannot perform dressing!");
143 return StatusCode::FAILURE;
146 listOfParticlesToDress.push_back(pItr);
150 bare_part.SetPtEtaPhiM(pt_visAcc(*pItr), eta_visAcc(*pItr), phi_visAcc(*pItr), mvisAcc(*pItr));
151 listOfDressedP4.push_back(bare_part);
155 listOfBareP4 = listOfDressedP4;
164 for (
const auto*
part : listOfParticlesToDress) {
165 listOfDressedP4.push_back(
part->p4());
169 std::vector<int> dressedParticlesNPhot(listOfParticlesToDress.size(), 0);
172 std::vector<const xAOD::TruthParticle*> photonsFSRList;
173 std::vector<int> photonPID{22};
175 photonPID, m_usePhotonsFromHadrons);
177 ATH_MSG_WARNING(
"Cannot construct the list of final state particles "<<m_truthClassKey.fullKey());
183 for (
const auto* phot : photonsFSRList) {
184 double dRmin = m_coneSize;
187 for (
size_t i = 0;
i < listOfParticlesToDress.size(); ++
i) {
188 if (!m_useLeptonsFromHadrons) {
189 if (!acc_origin.isAvailable(*listOfParticlesToDress[
i])) {
190 ATH_MSG_WARNING(
"MCTruthClassifier "<<m_truthClassKey.fullKey() <<
" not available, cannot apply notFromHadron veto!");
192 unsigned int result = acc_origin(*listOfParticlesToDress[
i]);
197 if (listOfParticlesToDress[
i]->
isTau()) {
198 bare_part = listOfBareP4[
i];
201 bare_part = listOfParticlesToDress[
i]->p4();
204 double dR = bare_part.DeltaR(phot->p4());
212 listOfDressedP4[
idx] += phot->p4();
213 dressedParticlesNPhot[
idx]++;
214 if (!m_decorationName.empty()){
221 for (
size_t i = 0;
i < listOfParticlesToDress.size(); ++
i) {
226 decorator_pt_vis(*
part) = dressedVec.Pt();
227 decorator_eta_vis(*
part) = dressedVec.Eta();
228 decorator_phi_vis(*
part) = dressedVec.Phi();
229 decorator_m_vis(*
part) = dressedVec.M();
232 decorator_e(*
part) = dressedVec.E();
233 decorator_pt(*
part) = dressedVec.Pt();
234 decorator_eta(*
part) = dressedVec.Eta();
235 decorator_phi(*
part) = dressedVec.Phi();
237 decorator_nphoton(*
part) = dressedParticlesNPhot[
i];
243 std::vector<fastjet::PseudoJet> sorted_jets;
244 std::vector<fastjet::PseudoJet> fj_particles;
245 for (
const auto*
part : listOfParticlesToDress) {
248 fj_particles.emplace_back(tauvis.Px(), tauvis.Py(), tauvis.Pz(), tauvis.E());
256 for (
const auto*
part : photonsFSRList) {
263 const fastjet::JetDefinition jet_def(
alg, m_coneSize);
265 sorted_jets = sorted_by_pt(cseq.inclusive_jets(0));
267 std::vector<int> photon_uniqueIDs(50);
268 photon_uniqueIDs.clear();
269 for (
const auto*
part : listOfParticlesToDress) {
272 auto pjItr=sorted_jets.begin();
273 while(!
found && pjItr!=sorted_jets.end()) {
274 std::vector<fastjet::PseudoJet> constituents = pjItr->constituents();
275 for(
const auto& constit : constituents) {
282 decorator_pt_vis(*
part) = pjItr->pt();
283 decorator_eta_vis(*
part) = pjItr->pseudorapidity();
284 decorator_phi_vis(*
part) = pjItr->phi_std();
285 decorator_m_vis(*
part) = pjItr->m();
288 decorator_e(*
part) = pjItr->e();
289 decorator_pt(*
part) = pjItr->pt();
290 decorator_eta(*
part) = pjItr->pseudorapidity();
291 decorator_phi(*
part) = pjItr->phi_std();
298 for(
const auto& constit : constituents) {
299 photon_uniqueIDs.push_back(constit.user_index());
306 decorator_pt_vis(*
part) = 0.;
307 decorator_eta_vis(*
part) = 0.;
308 decorator_phi_vis(*
part) = 0.;
309 decorator_m_vis(*
part) = 0.;
312 decorator_e(*
part) = 0;
313 decorator_pt(*
part) = 0;
314 decorator_eta(*
part) = 0;
315 decorator_phi(*
part) = 0;
321 if (!m_decorationName.empty()){
323 for (
const auto* phot : photonsFSRList ) {
332 return StatusCode::SUCCESS;