From f94e0f9f215edde703924865abe065d712bcb64a Mon Sep 17 00:00:00 2001
From: Christian Ohm <christian.ohm@cern.ch>
Date: Tue, 8 Jun 2021 00:25:05 +0200
Subject: [PATCH] Adding extra truth info to SiHits, optional of course

---
 .../Tools/HitAnalysis/src/SiHitAnalysis.cxx   | 43 ++++++++++++++++++-
 .../Tools/HitAnalysis/src/SiHitAnalysis.h     | 12 ++++--
 2 files changed, 51 insertions(+), 4 deletions(-)

diff --git a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx
index 317a0ec754d..058e006641e 100644
--- a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx
+++ b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.cxx
@@ -7,6 +7,9 @@ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
 // Section of includes for Pixel and SCT tests
 #include "InDetSimEvent/SiHitCollection.h"
 #include "GeoAdaptors/GeoSiHit.h"
+#include "GeneratorObjects/HepMcParticleLink.h"
+#include "HepMC/GenVertex.h"
+#include "HepMC/GenParticle.h"
 
 #include "TH1.h"
 #include "TH2.h"
@@ -42,9 +45,17 @@ SiHitAnalysis::SiHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator)
    , m_hits_eloss(0)
    , m_hits_step(0)
    , m_hits_barcode(0)
+   , m_hits_pdgId(0)
+   , m_hits_pT(0)
+   , m_hits_eta(0)
+   , m_hits_phi(0)
+   , m_hits_prodVtx_x(0)
+   , m_hits_prodVtx_y(0)
+   , m_hits_prodVtx_z(0)
      
    , m_collection("BCMHits")
    , m_expert("off")
+   , m_extraTruthBranches(0)
      
    , m_tree(0)
    , m_ntupleFileName("/SiHitAnalysis/")
@@ -59,6 +70,7 @@ SiHitAnalysis::SiHitAnalysis(const std::string& name, ISvcLocator* pSvcLocator)
   declareProperty("HistPath", m_path);
   declareProperty("isITK", m_isITK);
   declareProperty("isHGTD", m_isHGTD);
+  declareProperty("ExtraTruthBranches", m_extraTruthBranches = false);
 }
 
 StatusCode SiHitAnalysis::initialize() {
@@ -205,6 +217,15 @@ StatusCode SiHitAnalysis::initialize() {
     m_tree->Branch((detName+"_eloss").c_str(), &m_hits_eloss);
     m_tree->Branch((detName+"_step").c_str(), &m_hits_step);
     m_tree->Branch((detName+"_barcode").c_str(), &m_hits_barcode);
+    if (m_extraTruthBranches) {
+      m_tree->Branch((detName+"_pdgId").c_str(), &m_hits_pdgId);
+      m_tree->Branch((detName+"_pT").c_str(), &m_hits_pT);
+      m_tree->Branch((detName+"_eta").c_str(), &m_hits_eta);
+      m_tree->Branch((detName+"_phi").c_str(), &m_hits_phi);
+      m_tree->Branch((detName+"_prodVtx_x").c_str(), &m_hits_prodVtx_x);
+      m_tree->Branch((detName+"_prodVtx_y").c_str(), &m_hits_prodVtx_y);
+      m_tree->Branch((detName+"_prodVtx_z").c_str(), &m_hits_prodVtx_z);
+    }
   }
   else {
     ATH_MSG_ERROR("No tree found!");
@@ -225,6 +246,15 @@ StatusCode SiHitAnalysis::execute() {
   m_hits_eloss->clear();
   m_hits_step->clear();
   m_hits_barcode->clear();
+  if (m_extraTruthBranches) {
+    m_hits_pdgId->clear();
+    m_hits_pT->clear();
+    m_hits_eta->clear();
+    m_hits_phi->clear();
+    m_hits_prodVtx_x->clear();
+    m_hits_prodVtx_y->clear();
+    m_hits_prodVtx_z->clear();
+  }
   
   const DataHandle<SiHitCollection> p_collection;
   if (evtStore()->retrieve(p_collection, m_collection) == StatusCode::SUCCESS) {
@@ -261,7 +291,18 @@ StatusCode SiHitAnalysis::execute() {
       m_hits_eloss->push_back(i_hit->energyLoss());
       m_hits_time->push_back(i_hit->meanTime()); 
       m_hits_step->push_back(step_length);
-      m_hits_barcode->push_back(i_hit->particleLink().barcode());    
+      m_hits_barcode->push_back(i_hit->particleLink().barcode());
+      if (m_extraTruthBranches) {
+	auto tpl = i_hit->particleLink();
+	m_hits_pdgId->push_back(tpl->pdg_id());
+	m_hits_pT->push_back(tpl->momentum().perp());
+	m_hits_pT->push_back(tpl->momentum().eta());
+	m_hits_pT->push_back(tpl->momentum().phi());
+	m_hits_prodVtx_x->push_back(tpl->production_vertex()->position().x());
+	m_hits_prodVtx_y->push_back(tpl->production_vertex()->position().y());
+	m_hits_prodVtx_z->push_back(tpl->production_vertex()->position().z());
+      }
+
     } // End while hits
   } // End statuscode success upon retrieval of hits
 
diff --git a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h
index 4b72efb8f5a..c6bd2f33cb7 100755
--- a/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h
+++ b/Simulation/Tools/HitAnalysis/src/SiHitAnalysis.h
@@ -19,7 +19,6 @@ class TH1;
 class TH2;
 class TTree;
 
-
 class SiHitAnalysis : public AthAlgorithm {
 
  public:
@@ -55,9 +54,16 @@ class SiHitAnalysis : public AthAlgorithm {
    std::vector<float>* m_hits_eloss;
    std::vector<float>* m_hits_step;
    std::vector<float>* m_hits_barcode;
-
+   std::vector<float>* m_hits_pdgId;
+   std::vector<float>* m_hits_pT;
+   std::vector<float>* m_hits_eta;
+   std::vector<float>* m_hits_phi;
+   std::vector<float>* m_hits_prodVtx_x;
+   std::vector<float>* m_hits_prodVtx_y;
+   std::vector<float>* m_hits_prodVtx_z;
    std::string m_collection;
-   std::string m_expert;  
+   std::string m_expert;
+   bool m_extraTruthBranches;
    
    TTree* m_tree;
    std::string m_ntupleFileName; 
-- 
GitLab