Tool to calculate track truth. More...

#include <MuonTrackTruthTool.h>

Inheritance diagram for Muon::MuonTrackTruthTool:

Collaboration diagram for Muon::MuonTrackTruthTool:

Classes
struct	SortResultByMatchedHits

Public Member Functions
	~MuonTrackTruthTool ()=default
	destructor More...

StatusCode	initialize ()
	AlgTool initilize. More...

ResultVec	match (const TruthTree &truth_tree, const TrackCollection &tracks) const
	perform truth matching for a given set of tracks More...

SegmentResultVec	match (const TruthTree &truth_tree, const std::vector< const MuonSegment * > &segments) const
	perform truth matching for a given set of segments More...

MuonTrackTruth	getTruth (const TruthTree &truth_tree, const Trk::Track &track, bool restrictedTruth=false) const
	get track truth More...

MuonTrackTruth	getTruth (const TruthTree &truth_tree, const std::vector< const MuonSegment * > &segments, bool restrictedTruth=false) const
	get segment truth for a list of segments, the segments will be considered to belong to the same muon More...

MuonTrackTruth	getTruth (const TruthTree &truth_tree, const Muon::MuonSegment &segment) const
	get segment truth More...

MuonTrackTruth	getTruth (const TruthTree &truth_tree, const std::vector< const Trk::MeasurementBase * > &measurements, bool restrictedTruth=false) const
	get truth for a give set of hits. More...

const TruthTree	createTruthTree (const TrackRecordCollection truthTrackCol, const McEventCollection mcEventCollection, const std::vector< const MuonSimDataCollection * > &muonSimData, const CscSimDataCollection *cscSimDataMap) const
	create truth tree from sim data More...

HepMC::ConstGenParticlePtr	getMother (const TruthTrajectory &traj, const int barcodeIn) const
	Returns the mother particle of the particle with barcodeIn if it is found in the truth trajectory. More...

HepMC::ConstGenParticlePtr	getAncestor (const TruthTrajectory &traj, const int barcodeIn) const
	Returns the ancestor particle of the particle with barcodeIn if it is found in the truth trajectory. More...

HepMC::ConstGenParticlePtr	getInitial (const TruthTrajectory &traj, const int barcodeIn) const
	Returns the initial particle of the particle with barcodeIn if it is found in the truth trajectory. More...

unsigned int	getNumberOfScatters (const TruthTrajectory &traj, const int barcodeIn) const
	Returns the number of steps a particle took while maintaining its PDG ID. More...

Private Member Functions
MuonTrackTruth	getTruth (const std::vector< const Trk::MeasurementBase * > &measurements, const TruthTreeEntry &truthEntry, bool restrictedTruth) const

void	addSimDataToTree (TruthTree &truth_tree, std::map< int, int > &barcode_map, const MuonSimDataCollection *simDataCol) const

void	addCscSimDataToTree (TruthTree &truth_tree, std::map< int, int > &barcode_map, const CscSimDataCollection *simDataCol) const

void	addMdtTruth (MuonTechnologyTruth &trackTruth, const Identifier &id, const Trk::MeasurementBase &meas, const MuonSimDataCollection &simCol) const

void	addClusterTruth (MuonTechnologyTruth &trackTruth, const Identifier &id, const Trk::MeasurementBase &meas, const MuonSimDataCollection &simCol) const

void	addClusterTruth (MuonTechnologyTruth &trackTruth, const Identifier &id, const Trk::MeasurementBase &meas, const CscSimDataCollection &simCol) const

void	addMissedHits (MuonTechnologyTruth &truth, const std::set< Identifier > &ids, const std::set< Identifier > &chids, const MuonSimDataCollection &simCol, bool restrictedTruth) const

void	addMissedHits (MuonTechnologyTruth &truth, const std::set< Identifier > &ids, const std::set< Identifier > &chids, const CscSimDataCollection &simCol, bool restrictedTruth) const

int	manipulateBarCode (int barcode) const

bool	selectPdg (int pdg) const

const std::pair< HepMC::ConstGenParticlePtr, unsigned int >	getInitialPair (const TruthTrajectory &traj, const int barcodeIn) const
	Returns the initial particle of the particle with barcodeIn if it is found in the truth trajectory. More...

Private Attributes
ServiceHandle< Muon::IMuonIdHelperSvc >	m_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}

PublicToolHandle< Muon::MuonEDMPrinterTool >	m_printer {this, "Printer", "Muon::MuonEDMPrinterTool/MuonEDMPrinterTool"}

ToolHandle< Trk::ITruthTrajectoryBuilder >	m_truthTrajectoryBuilder

Gaudi::Property< bool >	m_manipulateBarCode {this, "ManipulateBarCode", false}

Gaudi::Property< bool >	m_doSummary {this, "DoSummary", false}

Gaudi::Property< bool >	m_matchAllParticles {this, "MatchAllParticles", true}

Gaudi::Property< unsigned int >	m_minHits {this, "MinHits", 4}

Gaudi::Property< std::vector< int > >	m_pdgsToBeConsidered {this, "ConsideredPDGs", {}}

std::set< int >	m_selectedPdgs

Detailed Description

Tool to calculate track truth.

Definition at line 43 of file MuonTrackTruthTool.h.

Constructor & Destructor Documentation

◆ ~MuonTrackTruthTool()

Muon::MuonTrackTruthTool::~MuonTrackTruthTool ( )

default

destructor

Member Function Documentation

◆ addClusterTruth() [1/2]

void Muon::MuonTrackTruthTool::addClusterTruth	(	MuonTechnologyTruth &	trackTruth,
		const Identifier &	id,
		const Trk::MeasurementBase &	meas,
		const CscSimDataCollection &	simCol
	)		const

private

Definition at line 616 of file MuonTrackTruthTool.cxx.

                                                                                        {
         Identifier layid = m_idHelperSvc->layerId(id);
         Identifier chid = m_idHelperSvc->chamberId(id);
         // find SimData corresponding to identifier
         CscSimDataCollection::const_iterator it = simCol.find(id);
         if (it == simCol.end()) {
             truth.wrongHits.insert(layid);
             truth.wrongChambers.insert(chid);
             return;
         }
  
         std::vector<CscSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
         std::vector<CscSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
         for (; dit != dit_end; ++dit) {
             truth.matchedHits.insert(layid);
             truth.matchedChambers.insert(chid);
         }
     }

◆ addClusterTruth() [2/2]

void Muon::MuonTrackTruthTool::addClusterTruth	(	MuonTechnologyTruth &	trackTruth,
		const Identifier &	id,
		const Trk::MeasurementBase &	meas,
		const MuonSimDataCollection &	simCol
	)		const

private

Definition at line 560 of file MuonTrackTruthTool.cxx.

                                                                                         {
         Identifier layid = m_idHelperSvc->layerId(id);
         Identifier chid = m_idHelperSvc->chamberId(id);
  
         MuonSimDataCollection::const_iterator it = simCol.end();
  
         bool goodCluster = false;
         const CompetingMuonClustersOnTrack* crot = dynamic_cast<const CompetingMuonClustersOnTrack*>(&meas);
         if (crot) {
             for (unsigned int i = 0; i < crot->numberOfContainedROTs(); ++i) {
                 const MuonClusterOnTrack* cluster = &crot->rioOnTrack(i);
                 if (!cluster) continue;
                 // find SimData corresponding to identifier
                 it = simCol.find(cluster->identify());
                 if (it != simCol.end()) {
                     goodCluster = true;
                     break;
                 }
             }
         } else {
             // Find SimData corresponding to identifier
             const Trk::RIO_OnTrack* rot = nullptr;
             Trk::RoT_Extractor::extract(rot, &meas);
             const Trk::PrepRawData* prd = rot->prepRawData();
             if (prd) {
                 // check if an identifier from the list of RDOs is matched to that in the SDO collection
                 const std::vector<Identifier> &rdoList = prd->rdoList();
                 std::vector<Identifier>::const_iterator rit = rdoList.begin();
                 std::vector<Identifier>::const_iterator rit_end = rdoList.end();
                 for (; rit != rit_end; ++rit) {
                     it = simCol.find(*rit);
                     if (it != simCol.end()) {
                         goodCluster = true;
                         break;
                     }
                 }
             } else {
                 it = simCol.find(id);
                 if (it != simCol.end()) goodCluster = true;
             }
         }
  
         if (!goodCluster || it == simCol.end()) {
             truth.wrongHits.insert(id);
             return;
         }
  
         std::vector<MuonSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
         std::vector<MuonSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
         for (; dit != dit_end; ++dit) {
             truth.matchedHits.insert(layid);
             truth.matchedChambers.insert(chid);
         }
     }

◆ addCscSimDataToTree()

void Muon::MuonTrackTruthTool::addCscSimDataToTree	(	TruthTree &	truth_tree,
		std::map< int, int > &	barcode_map,
		const CscSimDataCollection *	simDataCol
	)		const

private

Definition at line 319 of file MuonTrackTruthTool.cxx.

                                                                                                {
         // loop over sim collection and check whether identifiers are on track
         CscSimDataCollection::const_iterator it = simDataCol->begin();
         CscSimDataCollection::const_iterator it_end = simDataCol->end();
         for (; it != it_end; ++it) {
             Identifier id = it->first;
  
             // loop over deposits
             std::vector<CscSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
             std::vector<CscSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
             for (; dit != dit_end; ++dit) {
                 int barcodeIn = HepMC::uniqueID(dit->first);
                 std::map<int, int>::const_iterator bit = barcode_map.find(barcodeIn);
                 if (bit == barcode_map.end()) {
                     ATH_MSG_VERBOSE(" discarding "
                                     << "  " << m_idHelperSvc->toString(id) << "   barcode " << barcodeIn);
                     continue;
                 }
                 // replace barcode with barcode from map
                 int barcode = bit->second;
  
                 TruthTreeIt eit = truth_tree.find(barcode);
                 if (eit == truth_tree.end()) {
                     ATH_MSG_VERBOSE(" discarding "
                                     << "  " << m_idHelperSvc->toString(id) << "   barcode " << barcode);
                     continue;
                 }
  
                 if (msgLvl(MSG::VERBOSE)) {
                     ATH_MSG_VERBOSE(" adding hit " << m_idHelperSvc->toString(id) << "   barcode " << barcode);
                     if (barcode != barcodeIn) ATH_MSG_VERBOSE(" hit barcode " << barcodeIn);
                 }
                 eit->second.cscHits.insert(*it);
             }
         }
     }

◆ addMdtTruth()

void Muon::MuonTrackTruthTool::addMdtTruth	(	MuonTechnologyTruth &	trackTruth,
		const Identifier &	id,
		const Trk::MeasurementBase &	meas,
		const MuonSimDataCollection &	simCol
	)		const

private

Definition at line 525 of file MuonTrackTruthTool.cxx.

                                                                                     {
         const MdtDriftCircleOnTrack* mdt = dynamic_cast<const MdtDriftCircleOnTrack*>(&meas);
         if (!mdt) {
             ATH_MSG_WARNING(" dynamic_cast to MdtDriftCircleOnTrack failed for measurement with id " << m_idHelperSvc->toString(id));
             return;
         }
  
         Identifier chid = m_idHelperSvc->chamberId(id);
  
         // find SimData corresponding to identifier
         MuonSimDataCollection::const_iterator it = simCol.find(id);
         if (it == simCol.end()) {
             truth.wrongHits.insert(id);
             truth.wrongChambers.insert(chid);
             return;
         }
  
         std::vector<MuonSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
         std::vector<MuonSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
         for (; dit != dit_end; ++dit) {
             double radius = dit->second.firstEntry();
  
             // check whether SDO and hit have same sign, don't check sign if drift radius < 0.3 mm as it will be badly determined
             double checkSign = fabs(mdt->driftRadius()) > 0.3 ? radius * mdt->driftRadius() : 1;
             if (checkSign >= 0) {
                 truth.matchedHits.insert(id);
                 truth.matchedChambers.insert(chid);
             } else {
                 truth.wrongHits.insert(id);
                 truth.wrongChambers.insert(chid);
             }
         }
     }

◆ addMissedHits() [1/2]

void Muon::MuonTrackTruthTool::addMissedHits	(	MuonTechnologyTruth &	truth,
		const std::set< Identifier > &	ids,
		const std::set< Identifier > &	chids,
		const CscSimDataCollection &	simCol,
		bool	restrictedTruth
	)		const

private

Definition at line 498 of file MuonTrackTruthTool.cxx.

                                                                                                            {
         // loop over sim collection and check whether identifiers are on track
         CscSimDataCollection::const_iterator it = simCol.begin();
         CscSimDataCollection::const_iterator it_end = simCol.end();
         for (; it != it_end; ++it) {
             Identifier id = m_idHelperSvc->layerId(it->first);
  
             int isOnTrack = ids.count(id);
             if (isOnTrack) continue;
  
             // if restricted truth mode, skip if chamber has not hits
             Identifier chid = m_idHelperSvc->chamberId(id);
             bool chamberHasHits = chids.count(chid);
             if (restrictedTruth && !chamberHasHits) continue;
  
             // loop over deposits
             std::vector<CscSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
             std::vector<CscSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
             for (; dit != dit_end; ++dit) {
                 truth.missedHits.insert(id);
                 // only fill a missed chamber if the chamber had no hits
                 if (!chamberHasHits) truth.missedChambers.insert(chid);
             }
         }
     }

◆ addMissedHits() [2/2]

void Muon::MuonTrackTruthTool::addMissedHits	(	MuonTechnologyTruth &	truth,
		const std::set< Identifier > &	ids,
		const std::set< Identifier > &	chids,
		const MuonSimDataCollection &	simCol,
		bool	restrictedTruth
	)		const

private

Definition at line 469 of file MuonTrackTruthTool.cxx.

                                                                                                             {
         // loop over sim collection and check whether identifiers are on track
         MuonSimDataCollection::const_iterator it = simCol.begin();
         MuonSimDataCollection::const_iterator it_end = simCol.end();
         for (; it != it_end; ++it) {
             Identifier id = it->first;
             // for trigger chambers use layer id
             if (m_idHelperSvc->isTrigger(id) || m_idHelperSvc->isCsc(id)) id = m_idHelperSvc->layerId(id);
  
             int isOnTrack = ids.count(id);
             if (isOnTrack) continue;
  
             // if restricted truth mode, skip if chamber has not hits
             Identifier chid = m_idHelperSvc->chamberId(id);
             bool chamberHasHits = chids.count(chid);
             if (restrictedTruth && !chamberHasHits) continue;
  
             // loop over deposits
             std::vector<MuonSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
             std::vector<MuonSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
             for (; dit != dit_end; ++dit) {
                 truth.missedHits.insert(id);
                 // only fill a missed chamber if the chamber had no hits
                 if (!chamberHasHits) truth.missedChambers.insert(chid);
             }
         }
     }

◆ addSimDataToTree()

void Muon::MuonTrackTruthTool::addSimDataToTree	(	TruthTree &	truth_tree,
		std::map< int, int > &	barcode_map,
		const MuonSimDataCollection *	simDataCol
	)		const

private

Definition at line 266 of file MuonTrackTruthTool.cxx.

                                                                                              {
         // loop over sim collection and check whether identifiers are on track
         MuonSimDataCollection::const_iterator it = simDataCol->begin();
         MuonSimDataCollection::const_iterator it_end = simDataCol->end();
         for (; it != it_end; ++it) {
             Identifier id = it->first;
  
             // loop over deposits
             std::vector<MuonSimData::Deposit>::const_iterator dit = it->second.getdeposits().begin();
             std::vector<MuonSimData::Deposit>::const_iterator dit_end = it->second.getdeposits().end();
             for (; dit != dit_end; ++dit) {
                 int barcodeIn = HepMC::barcode(dit->first); // FIXME barcode-based
                 std::map<int, int>::const_iterator bit = barcode_map.find(barcodeIn);
                 if (bit == barcode_map.end()) {
                     ATH_MSG_VERBOSE(" discarding "
                                     << "  " << m_idHelperSvc->toString(id) << "   barcode " << barcodeIn);
                     continue;
                 }
                 // replace barcode with barcode from map
                 int barcode = bit->second;
  
                 TruthTreeIt eit = truth_tree.find(barcode);
                 if (eit == truth_tree.end()) {
                     ATH_MSG_VERBOSE(" discarding "
                                     << "  " << m_idHelperSvc->toString(id) << "   barcode " << barcode);
                     continue;
                 }
  
                 if (m_idHelperSvc->isMdt(id)) {
                     eit->second.mdtHits.insert(*it);
                 } else if (m_idHelperSvc->isRpc(id)) {
                     if (m_idHelperSvc->stationIndex(id) == StIndex::BO && m_idHelperSvc->rpcIdHelper().doubletR(id) == 2) {
                         ATH_MSG_VERBOSE(" Discarding non existing RPC hit " << m_idHelperSvc->toString(id));
                         continue;
                     }
  
                     eit->second.rpcHits.insert(*it);
                 } else if (m_idHelperSvc->isTgc(id)) {
                     eit->second.tgcHits.insert(*it);
                 } else if (m_idHelperSvc->issTgc(id)) {
                     eit->second.stgcHits.insert(*it);
                 } else if (m_idHelperSvc->isMM(id)) {
                     eit->second.mmHits.insert(*it);
                 }
                 if (msgLvl(MSG::VERBOSE)) {
                     ATH_MSG_VERBOSE(" adding hit " << m_idHelperSvc->toString(id) << "   barcode " << barcode);
                     if (barcode != barcodeIn) ATH_MSG_VERBOSE(" hit barcode " << barcodeIn);
                 }
             }
         }
     }

◆ createTruthTree()

const MuonTrackTruthTool::TruthTree Muon::MuonTrackTruthTool::createTruthTree	(	const TrackRecordCollection *	truthTrackCol,
		const McEventCollection *	mcEventCollection,
		const std::vector< const MuonSimDataCollection * > &	muonSimData,
		const CscSimDataCollection *	cscSimDataMap
	)		const

create truth tree from sim data

Definition at line 108 of file MuonTrackTruthTool.cxx.

                                                                                                                              {
         std::map<int, int> barcode_map;
         MuonTrackTruthTool::TruthTree truth_tree;
         if (truthTrackCol->empty()) {
             ATH_MSG_WARNING(" TrackRecordCollection is empty ");
             return truth_tree;
         }
  
         const HepMC::GenEvent* genEvent = nullptr;
         if (!mcEventCollection->empty()) {
             ATH_MSG_VERBOSE("McEventCollection size " << mcEventCollection->size());
             if (mcEventCollection->size() == 1) genEvent = mcEventCollection->front();
         }
  
         ATH_MSG_VERBOSE(" creating truth tree from track record " << truthTrackCol->size());
  
         TrackRecordConstIterator tr_it = truthTrackCol->begin();
         TrackRecordConstIterator tr_it_end = truthTrackCol->end();
         for (; tr_it != tr_it_end; ++tr_it) {
             int PDGCode((*tr_it).GetPDGCode());
             int barcode = HepMC::uniqueID(*tr_it);
             if (!m_matchAllParticles && !selectPdg(PDGCode)) {
                 ATH_MSG_VERBOSE(" discarding truth track: pdg " << PDGCode << "  barcode " << barcode);
                 continue;
             }
  
             // check whether barcode is already in, skip if that is the case
             if (barcode_map.count(barcode)) {
                 ATH_MSG_VERBOSE(" barcode " << barcode << " already in map, final state barcode " << barcode_map[barcode]);
                 continue;
             }
             ATH_MSG_VERBOSE(" found new particle with pdgid " << PDGCode << " in truth record, barcode " << barcode);
  
             std::unique_ptr<TruthTrajectory> truthTrajectory;
             // associate the muon truth with the gen event info
             if (genEvent) {
                 HepMC::ConstGenParticlePtr genParticle =
 #ifdef HEPMC3
                  genEvent->particles().at(HepMC::uniqueID(*tr_it)-1); // FIXME implement HepMC::id_to_particle/vertex explicitly
 #else
                  genEvent->barcode_to_particle(HepMC::uniqueID(*tr_it));
 #endif
                 if (genParticle) {
                     truthTrajectory = std::make_unique<TruthTrajectory>();
                     m_truthTrajectoryBuilder->buildTruthTrajectory(truthTrajectory.get(), genParticle);
                     if (!truthTrajectory->empty()) {
                         // always use barcode of the 'final' particle in chain in map
                         barcode = HepMC::uniqueID(truthTrajectory->front());
  
                         if (msgLvl(MSG::VERBOSE)) {
                             auto particle = truthTrajectory->front().cptr();
                             ATH_MSG_VERBOSE(" found GenParticle: size " << truthTrajectory->size() << " fs barcode " << barcode << "particle "<< particle);
                             if (particle->production_vertex()) {
                                 ATH_MSG_VERBOSE(" vertex: r  " << particle->production_vertex()->position().perp() << " z "
                                                                << particle->production_vertex()->position().z());
                             }
                         }
  
                         // now collect all barcodes beloning to this TruthTrajectory
                         std::vector<HepMcParticleLink>::const_iterator pit = truthTrajectory->begin();
                         std::vector<HepMcParticleLink>::const_iterator pit_end = truthTrajectory->end();
                         for (; pit != pit_end; ++pit) {
                             int code = HepMC::uniqueID(*pit);
  
                             if (msgLvl(MSG::VERBOSE) && code != barcode) {
                                 auto particle = (*pit).cptr();
                                 ATH_MSG_VERBOSE("  secondary barcode: " << code << "particle "<< particle);
                                 if (particle->production_vertex())
                                     ATH_MSG_VERBOSE(" vertex: r  " <<particle->production_vertex()->position().perp() << " z "
                                                                    << particle->production_vertex()->position().z());
                                 // sanity check
                                 if (barcode_map.count(code)) ATH_MSG_VERBOSE("  pre-existing barcode " << code);
                             }
  
                             // enter barcode
                             barcode_map[code] = barcode;
                         }
                     } else {
                         ATH_MSG_WARNING(" empty truth trajectory " << barcode);
                     }
                 }
             } else {
                 // add one to one relation
                 barcode_map[barcode] = barcode;
             }
  
             if (truth_tree.count(barcode)) {
                 ATH_MSG_WARNING(" found muon barcode twice in truth record: " << barcode);
                 continue;
             }
  
             TruthTreeEntry& entry = truth_tree[barcode];
             entry.truthTrack = &(*tr_it);
             // entry.truthTrajectory = truthTrajectory.get();
             entry.truthTrajectory = std::move(truthTrajectory);
             // m_truthTrajectoriesToBeDeleted.push_back(std::move(truthTrajectory));
         }
  
         // add sim data collections
         for (const MuonSimDataCollection* simDataMap : muonSimData) { addSimDataToTree(truth_tree, barcode_map, simDataMap); }
         if (cscSimDataMap) { addCscSimDataToTree(truth_tree, barcode_map, cscSimDataMap); }
  
         unsigned int ngood(0);
         std::vector<int> badBarcodes;
         // erase entries with too few hits or no track record
         TruthTreeIt it = truth_tree.begin();
         for (; it != truth_tree.end(); ++it) {
             bool erase = false;
             unsigned int nhits = it->second.mdtHits.size() + it->second.rpcHits.size() + it->second.tgcHits.size() +
                                  it->second.cscHits.size() + it->second.stgcHits.size() + it->second.mmHits.size();
             if (!it->second.truthTrack) erase = true;
             if (nhits < m_minHits) erase = true;
  
             if (erase) {
                 ATH_MSG_VERBOSE(" Erasing entry: barcode " << HepMC::uniqueID(it->second.truthTrack) << " hits " << nhits);
                 badBarcodes.push_back(it->first);
             } else {
                 ++ngood;
                 ATH_MSG_VERBOSE(" Keeping entry: barcode " << HepMC::uniqueID(it->second.truthTrack) << " hits " << nhits);
             }
         }
  
         std::vector<int>::iterator badIt = badBarcodes.begin();
         std::vector<int>::iterator badIt_end = badBarcodes.end();
         for (; badIt != badIt_end; ++badIt) truth_tree.erase(*badIt);
  
         if (ngood != truth_tree.size()) {
             ATH_MSG_WARNING(" Problem cleaning map: size " << truth_tree.size() << " accepted entries " << ngood);
         }
  
         if (m_doSummary || msgLvl(MSG::DEBUG)) {
             ATH_MSG_INFO(" summarizing truth tree: number of particles " << truth_tree.size());
             TruthTreeIt it = truth_tree.begin();
             TruthTreeIt it_end = truth_tree.end();
             for (; it != it_end; ++it) {
                 if (!it->second.truthTrack)
                     ATH_MSG_INFO(" no TrackRecord ");
                 else {
                     ATH_MSG_INFO(" PDG " << it->second.truthTrack->GetPDGCode() << " barcode " << HepMC::uniqueID(it->second.truthTrack));
                 }
                 if (!it->second.mdtHits.empty()) ATH_MSG_INFO(" mdt  " << it->second.mdtHits.size());
                 if (!it->second.rpcHits.empty()) ATH_MSG_INFO(" rpc  " << it->second.rpcHits.size());
                 if (!it->second.tgcHits.empty()) ATH_MSG_INFO(" tgc  " << it->second.tgcHits.size());
                 if (!it->second.cscHits.empty()) ATH_MSG_INFO(" csc  " << it->second.cscHits.size());
                 if (!it->second.stgcHits.empty()) ATH_MSG_INFO(" stgc " << it->second.stgcHits.size());
                 if (!it->second.mmHits.empty()) ATH_MSG_INFO(" mm   " << it->second.mmHits.size());
                 if (it->second.mdtHits.empty() && it->second.rpcHits.empty() && it->second.tgcHits.empty() && it->second.cscHits.empty() &&
                     it->second.stgcHits.empty() && it->second.mmHits.empty())
                     ATH_MSG_INFO(" no hits ");
             }
         }
  
         return truth_tree;
     }

◆ getAncestor()

HepMC::ConstGenParticlePtr Muon::MuonTrackTruthTool::getAncestor	(	const TruthTrajectory &	traj,
		const int	barcodeIn
	)		const

Returns the ancestor particle of the particle with barcodeIn if it is found in the truth trajectory.

Ancestor here means the last particle at generator level that has a status code different from final state, e.g. Z

Definition at line 658 of file MuonTrackTruthTool.cxx.

                                                                                                                    {
         bool foundBC = false;
         for (const auto& pit : traj) {
             if (!pit) continue;
             if (HepMC::uniqueID(pit) == barcodeIn || foundBC) {
                 foundBC = true;
 #ifdef HEPMC3
                 auto particle = pit.scptr();
 #else
                 auto particle = pit.cptr();
 #endif
                 if (!MC::isStable(particle)) {  // first non final state particle
                   return particle;
                 }
             }
         }
         return nullptr;
     }

◆ getInitial()

HepMC::ConstGenParticlePtr Muon::MuonTrackTruthTool::getInitial	(	const TruthTrajectory &	traj,
		const int	barcodeIn
	)		const

Returns the initial particle of the particle with barcodeIn if it is found in the truth trajectory.

For example a mu undergoing a mubrem would create a second mu, in which case this method returns the mu prior to bremsstrahlung. This interface calls the method getInitialPair.

Definition at line 737 of file MuonTrackTruthTool.cxx.

                                                                                                                   {
         return getInitialPair(traj, barcodeIn).first;
     }

◆ getInitialPair()

const std::pair< HepMC::ConstGenParticlePtr, unsigned int > Muon::MuonTrackTruthTool::getInitialPair	(	const TruthTrajectory &	traj,
		const int	barcodeIn
	)		const

private

Returns the initial particle of the particle with barcodeIn if it is found in the truth trajectory.

For example a mu undergoing a mubrem would create a second mu, in which case this method returns the mu prior to bremsstrahlung. The number of such scatters is returned in the .second.

Definition at line 677 of file MuonTrackTruthTool.cxx.

                                                                                                                             {
         std::pair<HepMC::ConstGenParticlePtr, unsigned int> thePair(nullptr, 0);
         unsigned int scat = 0;
         ATH_MSG_DEBUG("getFirst() : size = " << traj.size());
         bool foundBC = false;
         int pdgFinal = 0;
         double ePrev = 0.;
         HepMC::ConstGenParticlePtr theFirst{nullptr};
         for (auto pit = traj.begin(); pit != traj.end(); ++pit) {
            if (HepMC::uniqueID(*pit) == barcodeIn || foundBC) {
               auto particle = (*pit).scptr();
 #ifdef HEPMC3
                 if (!foundBC) {
                     foundBC = true;
                     theFirst = particle;
                     pdgFinal = particle->pdg_id();
                 } else {
                     if (particle->pdg_id() == pdgFinal) {
                         const auto& pit_p = *pit;
                         if ((theFirst != pit_p.scptr()) && (particle->momentum().t() != ePrev))
                             ++scat;  // if the particle has not changed pdgid after the first step count as scatter. also avoid counting
                                      // pure interface changes as scatter
                     } else {         // the first time this particle appears
                         --pit;       // AV This is confusing
                         theFirst = (*pit).scptr();
                         break;
                     }
                 }
 #else
                 if (!foundBC) {
                     foundBC = true;
                     theFirst = (*pit).cptr();
                     pdgFinal = (*pit)->pdg_id();
                 } else {
                     if ((*pit)->pdg_id() == pdgFinal) {
                         auto pit_p = *pit;
                         if ((theFirst != pit_p.cptr()) && ((*pit).cptr()->momentum().t() != ePrev))
                             ++scat;  // if the particle has not changed pdgid after the first step count as scatter. also avoid counting
                                      // pure interface changes as scatter
                     } else {         // the first time this particle appears
                         --pit;
                         theFirst = (*pit).cptr();
                         break;
                     }
                 }
 #endif
                 ATH_MSG_DEBUG("getFirst() : pt = " << particle->momentum().perp() << " scat = " << scat);
                 ePrev = particle->momentum().t();  // prepare for comparing this entry with the next one
             }
         }
         // sanity check
         if (theFirst && theFirst->pdg_id() != pdgFinal) ATH_MSG_ERROR("Wrong pdgId association in getFirst()");
         ATH_MSG_DEBUG("Number of scatters = " << scat << " pdgId = " << pdgFinal);
  
         thePair.first = theFirst;
         thePair.second = scat;
         return thePair;
     }

◆ getMother()

HepMC::ConstGenParticlePtr Muon::MuonTrackTruthTool::getMother	(	const TruthTrajectory &	traj,
		const int	barcodeIn
	)		const

Returns the mother particle of the particle with barcodeIn if it is found in the truth trajectory.

It traces the decay chain until if finds the first particle that is different flavor from the starting one.

Definition at line 636 of file MuonTrackTruthTool.cxx.

                                                                                                                  {
         ATH_MSG_DEBUG("getMother() : size = " << traj.size());
         int pdgFinal = ((traj.empty()) ? -999 : traj.front().cptr()->pdg_id());
         bool foundBC = false;
         for (const auto& pit : traj) {
             if (!pit) continue;
             if (HepMC::uniqueID(pit) == barcodeIn || foundBC) {
                 foundBC = true;
                 ATH_MSG_DEBUG("getMother() : " << pit );
 #ifdef HEPMC3
                 auto particle = pit.scptr();
 #else
                 auto particle = pit.cptr();
 #endif
                 if (particle->pdg_id() != pdgFinal) {  // the first case a track had a different flavour
                     break;
                 }
             }
         }
         return nullptr;
     }

◆ getNumberOfScatters()

unsigned int Muon::MuonTrackTruthTool::getNumberOfScatters	(	const TruthTrajectory &	traj,
		const int	barcodeIn
	)		const

Returns the number of steps a particle took while maintaining its PDG ID.

This method calls getInitialPair for calculating this number.

Definition at line 741 of file MuonTrackTruthTool.cxx.

                                                                                                                {
         return (getInitialPair(traj, barcodeIn)).second;
     }

◆ getTruth() [1/5]

MuonTrackTruth Muon::MuonTrackTruthTool::getTruth	(	const std::vector< const Trk::MeasurementBase * > &	measurements,
		const TruthTreeEntry &	truthEntry,
		bool	restrictedTruth
	)		const

private

Definition at line 421 of file MuonTrackTruthTool.cxx.

                                                                                                               {
         MuonTrackTruth trackTruth;
         trackTruth.truthTrack = truthEntry.truthTrack;
         trackTruth.truthTrajectory = truthEntry.truthTrajectory;
  
         std::vector<const Trk::MeasurementBase*>::const_iterator mit = measurements.begin();
         std::vector<const Trk::MeasurementBase*>::const_iterator mit_end = measurements.end();
         for (; mit != mit_end; ++mit) {
             // check whether state is a measurement
             const Trk::MeasurementBase* meas = *mit;
             if (!meas) { continue; }
  
             const Trk::RIO_OnTrack* rot = nullptr;
             Trk::RoT_Extractor::extract(rot, meas);
             if (!rot) {
                 if (!dynamic_cast<const Trk::PseudoMeasurementOnTrack*>(meas)) ATH_MSG_WARNING(" Could not get rot from measurement ");
                 continue;
             }
             Identifier id = rot->identify();
             if (!id.is_valid() || !m_idHelperSvc->mdtIdHelper().is_muon(id)) continue;
  
             if (m_idHelperSvc->isMdt(id)) {
                 addMdtTruth(trackTruth.mdts, id, *meas, truthEntry.mdtHits);
             } else if (m_idHelperSvc->isCsc(id)) {
                 addClusterTruth(trackTruth.cscs, id, *meas, truthEntry.cscHits);
             } else if (m_idHelperSvc->isRpc(id)) {
                 addClusterTruth(trackTruth.rpcs, id, *meas, truthEntry.rpcHits);
             } else if (m_idHelperSvc->isTgc(id)) {
                 addClusterTruth(trackTruth.tgcs, id, *meas, truthEntry.tgcHits);
             } else if (m_idHelperSvc->issTgc(id)) {
                 addClusterTruth(trackTruth.stgcs, id, *meas, truthEntry.stgcHits);
             } else if (m_idHelperSvc->isMM(id)) {
                 addClusterTruth(trackTruth.mms, id, *meas, truthEntry.mmHits);
             }
         }
  
         addMissedHits(trackTruth.mdts, trackTruth.mdts.matchedHits, trackTruth.mdts.matchedChambers, truthEntry.mdtHits, restrictedTruth);
         addMissedHits(trackTruth.cscs, trackTruth.cscs.matchedHits, trackTruth.cscs.matchedChambers, truthEntry.cscHits, restrictedTruth);
         addMissedHits(trackTruth.rpcs, trackTruth.rpcs.matchedHits, trackTruth.rpcs.matchedChambers, truthEntry.rpcHits, restrictedTruth);
         addMissedHits(trackTruth.tgcs, trackTruth.tgcs.matchedHits, trackTruth.tgcs.matchedChambers, truthEntry.tgcHits, restrictedTruth);
         addMissedHits(trackTruth.stgcs, trackTruth.stgcs.matchedHits, trackTruth.stgcs.matchedChambers, truthEntry.stgcHits,
                       restrictedTruth);
         addMissedHits(trackTruth.mms, trackTruth.mms.matchedHits, trackTruth.mms.matchedChambers, truthEntry.mmHits, restrictedTruth);
  
         return trackTruth;
     }

◆ getTruth() [2/5]

MuonTrackTruth Muon::MuonTrackTruthTool::getTruth	(	const TruthTree &	truth_tree,
		const Muon::MuonSegment &	segment
	)		const

get segment truth

Definition at line 357 of file MuonTrackTruthTool.cxx.

                                                                                                                  {
         return getTruth(truth_tree, segment.containedMeasurements(), true);
     }

◆ getTruth() [3/5]

MuonTrackTruth Muon::MuonTrackTruthTool::getTruth	(	const TruthTree &	truth_tree,
		const std::vector< const MuonSegment * > &	segments,
		bool	restrictedTruth = `false`
	)		const

get segment truth for a list of segments, the segments will be considered to belong to the same muon

Definition at line 366 of file MuonTrackTruthTool.cxx.

                                                                             {
         std::set<Identifier> ids;
         std::vector<const Trk::MeasurementBase*> measurements;
         std::vector<const MuonSegment*>::const_iterator sit = segments.begin();
         std::vector<const MuonSegment*>::const_iterator sit_end = segments.end();
         for (; sit != sit_end; ++sit) {
             std::vector<const Trk::MeasurementBase*>::const_iterator mit = (*sit)->containedMeasurements().begin();
             std::vector<const Trk::MeasurementBase*>::const_iterator mit_end = (*sit)->containedMeasurements().end();
             for (; mit != mit_end; ++mit) {
                 const Trk::MeasurementBase* meas = *mit;
                 const Trk::RIO_OnTrack* rot = nullptr;
                 Trk::RoT_Extractor::extract(rot, meas);
                 if (!rot) {
                     if (!dynamic_cast<const Trk::PseudoMeasurementOnTrack*>(meas)) ATH_MSG_WARNING(" Could not get rot from measurement ");
                     continue;
                 }
                 Identifier id = rot->identify();
                 if (!id.is_valid() || !m_idHelperSvc->mdtIdHelper().is_muon(id)) continue;
                 if (ids.count(id)) continue;
                 measurements.push_back(meas);
                 ids.insert(id);
             }
         }
         return getTruth(truth_tree, measurements, restrictedTruth);
     }

◆ getTruth() [4/5]

MuonTrackTruth Muon::MuonTrackTruthTool::getTruth	(	const TruthTree &	truth_tree,
		const std::vector< const Trk::MeasurementBase * > &	measurements,
		bool	restrictedTruth = `false`
	)		const

get truth for a give set of hits.

If restrictedTruth is set to true only missed hits in chambers with hits will be counted.

Definition at line 393 of file MuonTrackTruthTool.cxx.

                                                                             {
         MuonTrackTruth bestMatch;
         bestMatch.truthTrack = nullptr;
         bestMatch.truthTrajectory = nullptr;
  
         unsigned int nmatchedHitsBest = 0;
         // loop over muons and match hits
         TruthTreeConstIt tit = truth_tree.begin();
         TruthTreeConstIt tit_end = truth_tree.end();
         for (; tit != tit_end; ++tit) {
             unsigned int nhits = tit->second.mdtHits.size() + tit->second.cscHits.size() + tit->second.rpcHits.size() +
                                  tit->second.tgcHits.size() + tit->second.stgcHits.size() + tit->second.mmHits.size();
             if (nhits == 0) continue;
  
             MuonTrackTruth trackTruth = getTruth(measurements, tit->second, restrictedTruth);
             unsigned int nmatchedHits = trackTruth.numberOfMatchedHits();
             ATH_MSG_DEBUG(" performed truth match for particle with barcode: " << tit->first << " overlap " << nmatchedHits << " fraction  "
                                                                                << (double)nmatchedHits / (double)nhits);
             if (nmatchedHits > 0 && nmatchedHits > nmatchedHitsBest) {
                 bestMatch = trackTruth;
                 nmatchedHitsBest = nmatchedHits;
             }
         }
  
         return bestMatch;
     }

◆ getTruth() [5/5]

MuonTrackTruth Muon::MuonTrackTruthTool::getTruth	(	const TruthTree &	truth_tree,
		const Trk::Track &	track,
		bool	restrictedTruth = `false`
	)		const

get track truth

Definition at line 361 of file MuonTrackTruthTool.cxx.

                                                                                                                               {
         if (track.measurementsOnTrack()) return getTruth(truth_tree, track.measurementsOnTrack()->stdcont(), restrictedTruth);
         return {};
     }

◆ initialize()

StatusCode Muon::MuonTrackTruthTool::initialize ( )

AlgTool initilize.

Definition at line 25 of file MuonTrackTruthTool.cxx.

                                               {
         ATH_CHECK(m_idHelperSvc.retrieve());
         ATH_CHECK(m_printer.retrieve());
         ATH_CHECK(m_truthTrajectoryBuilder.retrieve());
  
         // add muons
         if (m_pdgsToBeConsidered.value().empty()) {
             m_selectedPdgs.insert(13);
             m_selectedPdgs.insert(-13);
         } else {
             // add pdgs
             for (auto pdg : m_pdgsToBeConsidered.value()) { m_selectedPdgs.insert(pdg); }
             ATH_MSG_DEBUG(" PDG codes used for matching "<<m_selectedPdgs);
         }
         return StatusCode::SUCCESS;
     }

◆ manipulateBarCode()

int Muon::MuonTrackTruthTool::manipulateBarCode ( int barcode ) const

private

Definition at line 42 of file MuonTrackTruthTool.cxx.

                                                                { // FIXME This is obsolete now
         if (m_manipulateBarCode) return barcode % 10000;
         return barcode;
     }

◆ match() [1/2]

MuonTrackTruthTool::SegmentResultVec Muon::MuonTrackTruthTool::match	(	const TruthTree &	truth_tree,
		const std::vector< const MuonSegment * > &	segments
	)		const

perform truth matching for a given set of segments

Definition at line 89 of file MuonTrackTruthTool.cxx.

                                                                                                                         {
         SegmentResultVec result;
         result.reserve(segments.size());
  
         // loop over tracks and match each of them with truth
         std::vector<const MuonSegment*>::const_iterator sit = segments.begin();
         std::vector<const MuonSegment*>::const_iterator sit_end = segments.end();
         for (; sit != sit_end; ++sit) {
             // create truth association
             result.emplace_back(*sit, getTruth(truth_tree, **sit));
         }
  
         // sort result per muon and per number of matched hits
         std::stable_sort(result.begin(), result.end(), SortResultByMatchedHits());
  
         return result;
     }

◆ match() [2/2]

MuonTrackTruthTool::ResultVec Muon::MuonTrackTruthTool::match	(	const TruthTree &	truth_tree,
		const TrackCollection &	tracks
	)		const

perform truth matching for a given set of tracks

Definition at line 66 of file MuonTrackTruthTool.cxx.

                                                                                                                           {
         ResultVec result;
         result.reserve(tracks.size());
  
         // loop over tracks and match each of them with truth
         TrackCollection::const_iterator tit = tracks.begin();
         TrackCollection::const_iterator tit_end = tracks.end();
         for (; tit != tit_end; ++tit) {
             MuonTrackTruth match = getTruth(truth_tree, **tit);
             if (!match.truthTrack) continue;
  
             if (match.numberOfMatchedHits() == 0) continue;
  
             // create truth association
             result.emplace_back(*tit, match);
         }
  
         // sort result per muon and per number of matched hits
         std::stable_sort(result.begin(), result.end(), SortResultByMatchedHits());
  
         return result;
     }

◆ selectPdg()

bool Muon::MuonTrackTruthTool::selectPdg ( int pdg ) const

inlineprivate

Definition at line 127 of file MuonTrackTruthTool.h.

127 { return m_selectedPdgs.count(pdg); }

Member Data Documentation

◆ m_doSummary

Gaudi::Property<bool> Muon::MuonTrackTruthTool::m_doSummary {this, "DoSummary", false}

private

Definition at line 141 of file MuonTrackTruthTool.h.

◆ m_idHelperSvc

ServiceHandle<Muon::IMuonIdHelperSvc> Muon::MuonTrackTruthTool::m_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}

private

Definition at line 134 of file MuonTrackTruthTool.h.

◆ m_manipulateBarCode

Gaudi::Property<bool> Muon::MuonTrackTruthTool::m_manipulateBarCode {this, "ManipulateBarCode", false}

private

Definition at line 140 of file MuonTrackTruthTool.h.

◆ m_matchAllParticles

Gaudi::Property<bool> Muon::MuonTrackTruthTool::m_matchAllParticles {this, "MatchAllParticles", true}

private

Definition at line 142 of file MuonTrackTruthTool.h.

◆ m_minHits

Gaudi::Property<unsigned int> Muon::MuonTrackTruthTool::m_minHits {this, "MinHits", 4}

private

Definition at line 143 of file MuonTrackTruthTool.h.

◆ m_pdgsToBeConsidered

Gaudi::Property<std::vector<int> > Muon::MuonTrackTruthTool::m_pdgsToBeConsidered {this, "ConsideredPDGs", {}}

private

Definition at line 144 of file MuonTrackTruthTool.h.

◆ m_printer

PublicToolHandle<Muon::MuonEDMPrinterTool> Muon::MuonTrackTruthTool::m_printer {this, "Printer", "Muon::MuonEDMPrinterTool/MuonEDMPrinterTool"}

private

Definition at line 136 of file MuonTrackTruthTool.h.

◆ m_selectedPdgs

std::set<int> Muon::MuonTrackTruthTool::m_selectedPdgs

private

Definition at line 146 of file MuonTrackTruthTool.h.

◆ m_truthTrajectoryBuilder

ToolHandle<Trk::ITruthTrajectoryBuilder> Muon::MuonTrackTruthTool::m_truthTrajectoryBuilder

private

Initial value:

{

this, "TruthTrajectoryBuilder", "Muon::MuonDecayTruthTrajectoryBuilder/MuonDecayTruthTrajectoryBuilder"}

Definition at line 137 of file MuonTrackTruthTool.h.

The documentation for this class was generated from the following files:

Classes

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ ~MuonTrackTruthTool()

Member Function Documentation

◆ addClusterTruth() [1/2]

◆ addClusterTruth() [2/2]

◆ addCscSimDataToTree()

◆ addMdtTruth()

◆ addMissedHits() [1/2]

◆ addMissedHits() [2/2]

◆ addSimDataToTree()

◆ createTruthTree()

◆ getAncestor()

◆ getInitial()

◆ getInitialPair()

◆ getMother()

◆ getNumberOfScatters()

◆ getTruth() [1/5]

◆ getTruth() [2/5]

◆ getTruth() [3/5]

◆ getTruth() [4/5]

◆ getTruth() [5/5]

◆ initialize()

◆ manipulateBarCode()

◆ match() [1/2]

◆ match() [2/2]

◆ selectPdg()

Member Data Documentation

◆ m_doSummary

◆ m_idHelperSvc

◆ m_manipulateBarCode

◆ m_matchAllParticles

◆ m_minHits

◆ m_pdgsToBeConsidered

◆ m_printer

◆ m_selectedPdgs

◆ m_truthTrajectoryBuilder