Commit d741ed07 authored by Alexander Leopold's avatar Alexander Leopold Committed by Tadej Novak
Browse files

HGTD_Cluster truth classification tool

parent a3c65c91
/**
* Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration.
*
* @file HGTD_RecToolInterfaces/IHGTD_ClusterTruthTool.h
* @author Alexander Leopold <alexander.leopold@cern.ch>
* @date August, 2021
*
* @brief To assess the performance of track hit matching algorithms, it is
* necessary to access the truth origin of hits in HGTD that were associated to
* tracks. It uses both RDO and SDO collections to trace back the truth info.
*/
#ifndef IHGTD_CLUSTERTRUTHTOOL_H
#define IHGTD_CLUSTERTRUTHTOOL_H
#include "GaudiKernel/IAlgTool.h"
namespace HepMC {
class GenEvent;
}
class InDetSimDataCollection;
// namespace xAOD {
// class TruthParticle;
// }
// FIXME ugly, but this works, the above does not. fails with:
// error: conflicting declaration 'typedef class xAOD::TruthParticle_v1
// xAOD::TruthParticle'
// typedef TruthParticle_v1 TruthParticle;
namespace xAOD {
class TruthParticle_v1;
typedef TruthParticle_v1 TruthParticle;
} // namespace xAOD
namespace HGTD {
class HGTD_Cluster;
enum ClusterTruthOrigin {
UNIDENTIFIED = 0,
TRUTH_PARTICLE = 1, // originates from the tested truth particle
HARD_SCATTER = 2, // originates from a HS particle, but unclear which
PILEUP = 3, // originates from a pileup interaction
SECONDARY = 4 // originates from some secondary interaction
};
struct ClusterTruthInfo {
ClusterTruthOrigin origin;
/**
* @brief Shadowing means that a deposit was left by the truth
* particle, but it was not the first deposit and is thus not used for the
* time measurement -> I will have an incorrect time
*/
bool is_shadowed;
/**
* @brief A cluster is considered to be merged if more than one particle
* deposited energy in a given pad.
*/
bool is_merged;
};
static const InterfaceID
IID_IHGTD_ClusterTruthTool("HGTD::IHGTD_ClusterTruthTool", 1, 0);
class IHGTD_ClusterTruthTool : virtual public IAlgTool {
public:
static const InterfaceID& interfaceID();
/**
* @brief Checks if the cluster originates from the tested truth particle or
* if it can be assigned to the hard scatter event. If not, the truth
* information can be used to determine which category (ClusterTruthOrigin)
* the cluster comes from.
*
* @param [in] cluster Hit in HGTD, for which we want the truth information.
* @param [in] tp Truth particle that is potentially matched to the cluster.
* @param [in] sim_data SDO collection storing the truth links of the charged
* diodes.
* @param [in] hard_scatter_evnt If given, a cluster can be categorised as
* coming from the HS event, even if the direct match with the truth particle
* fails.
*
* @return Struct combining the relevant truth information.
*/
virtual ClusterTruthInfo
classifyCluster(const HGTD::HGTD_Cluster* cluster,
const xAOD::TruthParticle* tp,
const InDetSimDataCollection* sim_data,
const HepMC::GenEvent* hard_scatter_evnt = nullptr) const = 0;
};
/** Inline methods **/
inline const InterfaceID& IHGTD_ClusterTruthTool::interfaceID() {
return IID_IHGTD_ClusterTruthTool;
}
} // namespace HGTD
#endif // IHGTD_CLUSTERTRUTHTOOL_H
################################################################################
# Package: HGTD_TruthTools
################################################################################
# Declare the package name:
atlas_subdir( HGTD_TruthTools )
find_package( ROOT COMPONENTS Core Tree MathCore Hist RIO pthread )
atlas_add_component( HGTD_TruthTools
src/*.cxx
src/components/*.cxx
LINK_LIBRARIES AthenaBaseComps GaudiKernel
HGTD_RecToolInterfaces GeneratorObjects ${ROOT_LIBRARIES}
HGTD_PrepRawData InDetSimData xAODTruth)
/**
* Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration.
*
* @file HGTD_TruthTools/ClusterTruthTool.h
*
* @author Alexander Leopold <alexander.leopold@cern.ch>
* @author Noemi Calace <noemi.calace@cern.ch>
*
* @date August, 2021
*
* @brief Uses the same categorisation as for the HGTD TDR results, taken from:
* https://gitlab.cern.ch/aleopold/memleakhunt/-/blob/master/TrkExtrapolation/TrkExTools/src/TrackTimingExtensionAlg.cxx#L901
*/
#ifndef HGTD_TRUTHTOOLS_CLUSTERTRUTHTOOL_H
#define HGTD_TRUTHTOOLS_CLUSTERTRUTHTOOL_H
#include "AthenaBaseComps/AthAlgTool.h"
#include "HGTD_RecToolInterfaces/IHGTD_ClusterTruthTool.h"
class HGTD_ID;
namespace HGTD {
class ClusterTruthTool : virtual public IHGTD_ClusterTruthTool,
public AthAlgTool {
public:
ClusterTruthTool(const std::string&, const std::string&, const IInterface*);
virtual ~ClusterTruthTool() {}
/**
* @brief The InDetSimDataCollection is a map, connecting Identifiers from the
* RDOs to InDetSimData objects, that keep a vector of HepMcParticleLink,
* float pairs, where the float keeps the energy of the given energy deposit.
* For HGTD, we store the time of the deposit instead, and can then select the
* first deposit, since only this one would contribute to the time read out by
* the ASIC.
*
* @param [in] cluster Hit in HGTD, for which we want the truth information.
* @param [in] tp Truth particle that is potentially matched to the cluster.
* @param [in] sim_data SDO collection storing the truth links of the charged
* diodes.
* @param [in] hard_scatter_evnt If given, a cluster can be categorised as
* coming from the HS event, even if the direct match with the truth particle
* fails.
*
* @return Struct combining the relevant truth information.
*/
virtual ClusterTruthInfo classifyCluster(
const HGTD::HGTD_Cluster* cluster, const xAOD::TruthParticle* tp,
const InDetSimDataCollection* sim_data,
const HepMC::GenEvent* hard_scatter_evnt = nullptr) const override;
};
} // namespace HGTD
#endif // HGTD_TRUTHTOOLS_CLUSTERTRUTHTOOL_H
/**
* Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration.
*
* @file HGTD_TruthTools/src/ClusterTruthTool.cxx
*
* @author Alexander Leopold <alexander.leopold@cern.ch>
* @author Noemi Calace <noemi.calace@cern.ch>
*
* @date August, 2021
*
* @brief
*/
#include "HGTD_TruthTools/ClusterTruthTool.h"
#include "GeneratorObjects/McEventCollection.h"
#include "HGTD_PrepRawData/HGTD_Cluster.h"
#include "InDetSimData/InDetSimDataCollection.h"
#include "xAODTruth/TruthParticle.h"
#include "TLorentzVector.h"
HGTD::ClusterTruthTool::ClusterTruthTool(const std::string& t,
const std::string& n,
const IInterface* p)
: AthAlgTool(t, n, p) {}
HGTD::ClusterTruthInfo HGTD::ClusterTruthTool::classifyCluster(
const HGTD::HGTD_Cluster* cluster, const xAOD::TruthParticle* tp,
const InDetSimDataCollection* sim_data,
const HepMC::GenEvent* hard_scatter_evnt) const {
const std::vector<Identifier>& rdo_id_list = cluster->rdoList();
// keep record of the cluster origins and if they are shadowed or not
std::vector<std::pair<HGTD::ClusterTruthOrigin, bool>> shadowed_origins;
for (const auto& rdo_id : rdo_id_list) {
auto pos = sim_data->find(rdo_id);
// the InDetSimData contains a std::pair<HepMcParticleLink, float>, where
// the second entry in the pair holds the time of the SiChargedDiode
if (pos == sim_data->end()) {
ATH_MSG_WARNING(
"[ClusterTruthTool::classifyCluster] ID not found in SDO map");
// FIXME I should probably continue here already? otherwise I get an
// "empty" entry in shadowed_origins
}
// collect deposits, sorted with first deposit at start of map bu default
std::map<float, ClusterTruthOrigin> sorted_deposits;
//////////////////////////////
// the following is taken from 20.20 as is
for (const auto& deposit : pos->second.getdeposits()) {
const HepMcParticleLink& particle_link = deposit.first;
int barcode = particle_link.barcode();
// barcodes of 0 or above 200k are secondaries generated by GEANT4
// 0 is also used for detector noise, delta rays and random energy
// deposits
if (barcode == 0 || barcode > 200000) {
sorted_deposits.emplace(deposit.second, ClusterTruthOrigin::SECONDARY);
// FIXME: shouldn't I "continue" here?
}
// check for identity with original particle
const HepMC::GenParticle* gen_part = particle_link.cptr();
if (gen_part) {
TLorentzVector l4(gen_part->momentum().px(), gen_part->momentum().py(),
gen_part->momentum().pz(), gen_part->momentum().e());
// if the barcode is identical and spacial matching passes, then this
// deposit came from the tested truth particle
if (barcode == tp->barcode() && tp->p4().DeltaR(l4) < 0.05) {
sorted_deposits.emplace(deposit.second,
ClusterTruthOrigin::TRUTH_PARTICLE);
// if given, the parent event can be checked
} else if (hard_scatter_evnt and
gen_part->parent_event() == hard_scatter_evnt) {
sorted_deposits.emplace(deposit.second,
ClusterTruthOrigin::HARD_SCATTER);
// otherwise, a particle that was generated but doesn't come from the
// hard scatter is considered to be originating from a pileup
// interaction
} else {
sorted_deposits.emplace(deposit.second, ClusterTruthOrigin::PILEUP);
}
} else {
// if there is no gen particle, we can guess based on the event index
if (particle_link.eventIndex() == 0) {
sorted_deposits.emplace(deposit.second,
ClusterTruthOrigin::HARD_SCATTER);
} else {
// If the gen was not kept, we assume it is cut away from truth record
// and is pileup
sorted_deposits.emplace(deposit.second, ClusterTruthOrigin::PILEUP);
}
}
} // END lOOP over the deposits
// check for shadowing, which means that a deposit was left by the truth
// particle, but it was not the first deposit and is thus not used for the
// time measurement -> I will have an incorrect time
bool is_shadowed = false;
ClusterTruthOrigin current_origin;
std::map<float, ClusterTruthOrigin>::iterator elem =
sorted_deposits.begin();
for (; elem != sorted_deposits.end(); ++elem) {
if (elem == sorted_deposits.begin()) {
current_origin = elem->second;
} else {
// if one of the later deposits originates from the truth particle, then
// the hit I want to find was shadowed by something else, and I will
// reconstruct an incorrect time
if (current_origin != ClusterTruthOrigin::TRUTH_PARTICLE &&
elem->second == ClusterTruthOrigin::TRUTH_PARTICLE) {
is_shadowed = true;
}
}
} // END LOOP over the time sorted deposits
shadowed_origins.emplace_back(current_origin, is_shadowed);
} // END LOOP over RDO identifiers
HGTD::ClusterTruthInfo result;
if (shadowed_origins.size() == 0) {
ATH_MSG_WARNING("did not manage to understand any RDOs...");
result.origin = ClusterTruthOrigin::UNIDENTIFIED;
result.is_shadowed = false;
// A cluster is considered to be merged if more than one particle deposited
// energy in a given pad.
result.is_merged = false;
} else {
result.is_merged = false;
result.origin = shadowed_origins.at(0).first;
result.is_shadowed = shadowed_origins.at(0).second;
for (size_t i = 1; i < shadowed_origins.size(); ++i) {
if (shadowed_origins.at(i).first != result.origin) {
result.is_merged = true;
} else {
result.is_shadowed |= shadowed_origins.at(i).second;
}
}
}
return result;
}
#include "GaudiKernel/DeclareFactoryEntries.h"
#include "HGTD_TruthTools/ClusterTruthTool.h"
using namespace HGTD;
DECLARE_TOOL_FACTORY(ClusterTruthTool)
DECLARE_FACTORY_ENTRIES( HGTD_TruthTools )
{
DECLARE_TOOL( ClusterTruthTool )
}
#include "GaudiKernel/DeclareFactoryEntries.h"
#include "GaudiKernel/LoadFactoryEntries.h"
LOAD_FACTORY_ENTRIES(HGTD_TruthTools)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment