ATLAS Offline Software
Loading...
Searching...
No Matches
LArHitEMap Class Reference

#include <LArHitEMap.h>

Collaboration diagram for LArHitEMap:

Public Types

typedef std::vector< const LArDigit * > DIGMAP

Public Member Functions

 LArHitEMap ()=delete
 LArHitEMap (const LArOnOffIdMapping *cabling, const CaloCell_ID *cellid, const CaloDetDescrManager *cddMgr, bool digit=false)
 ~LArHitEMap ()=default
bool AddEnergy (const IdentifierHash index, const float energy, const float time)
bool AddEnergy (const Identifier cellid, const float energy, const float time)
bool BuildWindows (const McEventCollection *mcCollptr, float deta, float dphi, float ptmin)
size_t GetNbCells (void) const
const LArHitListGetCell (const unsigned int index) const
const std::vector< std::pair< float, float > > & GetTimeE (const IdentifierHash index) const
bool AddDigit (const LArDigit *digit)
const LArDigitGetDigit (unsigned int index) const

Private Attributes

std::vector< LArHitListm_emap
DIGMAP m_digmap
const LArOnOffIdMappingm_cabling
const CaloCell_IDm_calocell_id
const CaloDetDescrManagerm_cddMgr

Detailed Description

Definition at line 21 of file LArHitEMap.h.

Member Typedef Documentation

◆ DIGMAP

typedef std::vector<const LArDigit*> LArHitEMap::DIGMAP

Definition at line 25 of file LArHitEMap.h.

Constructor & Destructor Documentation

◆ LArHitEMap() [1/2]

LArHitEMap::LArHitEMap ( )
delete

◆ LArHitEMap() [2/2]

LArHitEMap::LArHitEMap ( const LArOnOffIdMapping * cabling,
const CaloCell_ID * cellid,
const CaloDetDescrManager * cddMgr,
bool digit = false )

Definition at line 22 of file LArHitEMap.cxx.

22 :
23 m_cabling(cabling),
24 m_calocell_id(cellid),
25 m_cddMgr(cddMgr) {
26
27 //the last cell of the FCAL is the hash-max for LAr (ignore the Tile part)
28 IdentifierHash fcalCellMin, fcalCellMax;
29 cellid->calo_cell_hash_range(CaloCell_ID::LARFCAL,fcalCellMin,fcalCellMax);
30
31 //fill energy map up to fcal-hashmax(= lar-hashmax)
32 m_emap.resize(fcalCellMax);
33 if (digit) m_digmap.resize(fcalCellMax,nullptr);
34}
void calo_cell_hash_range(const Identifier id, IdentifierHash &caloCellMin, IdentifierHash &caloCellMax) const
to loop on 'global' cell hashes of one sub-calorimeter alone
std::vector< LArHitList > m_emap
Definition LArHitEMap.h:28
const CaloCell_ID * m_calocell_id
Definition LArHitEMap.h:31
const LArOnOffIdMapping * m_cabling
Definition LArHitEMap.h:30
const CaloDetDescrManager * m_cddMgr
Definition LArHitEMap.h:32
DIGMAP m_digmap
Definition LArHitEMap.h:29

◆ ~LArHitEMap()

LArHitEMap::~LArHitEMap ( )
default

Member Function Documentation

◆ AddDigit()

bool LArHitEMap::AddDigit ( const LArDigit * digit)

Definition at line 50 of file LArHitEMap.cxx.

50 {
51 const HWIdentifier ch_id = digit->channelID();
52 if (m_cabling->isOnlineConnected(ch_id)) {
53 Identifier cellid=m_cabling->cnvToIdentifier(ch_id);
54 IdentifierHash h=m_calocell_id->calo_cell_hash(cellid);
55
56 if (h>=m_digmap.size()) return false ;
58 return true;
59 }
60 else
61 return false;
62}

◆ AddEnergy() [1/2]

bool LArHitEMap::AddEnergy ( const Identifier cellid,
const float energy,
const float time )

Definition at line 45 of file LArHitEMap.cxx.

45 {
46 IdentifierHash idHash=m_calocell_id->calo_cell_hash(cellid);
47 return AddEnergy(idHash,energy,time);
48}
bool AddEnergy(const IdentifierHash index, const float energy, const float time)

◆ AddEnergy() [2/2]

bool LArHitEMap::AddEnergy ( const IdentifierHash index,
const float energy,
const float time )

Definition at line 38 of file LArHitEMap.cxx.

38 {
39 if(index >= m_emap.size()) return(false);
40 m_emap[index].AddHit(energy,time);
41 return true;
42}
str index
Definition DeMoScan.py:362

◆ BuildWindows()

bool LArHitEMap::BuildWindows ( const McEventCollection * mcCollptr,
float deta,
float dphi,
float ptmin )

Definition at line 64 of file LArHitEMap.cxx.

66{
67// get list of particles
68 std::vector<double> phiPart;
69 std::vector<double> etaPart;
70
71 etaPart.clear();
72 phiPart.clear();
73
74 if (!mcCollptr) {
75 return false;
76 }
77
79// std::cout << " start loop over particles " << std::endl;
80 for (itr = mcCollptr->begin(); itr!=mcCollptr->end(); ++itr) {
81 for (const auto& part: *(*itr))
82 {
83 //works only for photons(22) and electrons(11) primary particle (+pi0 in case not decayed by generator)
84 // with pt>5 GeV
85 // pickup "stable" particle from generator excluding G4 secondaries
86 if( ( MC::isPhoton(part) || MC::isElectron(part) || part->pdg_id()==111) && part->momentum().perp()> ptmin
88 {
89 etaPart.push_back(part->momentum().pseudoRapidity());
90 phiPart.push_back(part->momentum().phi());
91 }
92 }
93 }
94
95
96 if ( etaPart.empty()) return true;
97
98 for (unsigned int i=0; i < m_emap.size(); i++)
99 {
100 LArHitList& theLArHitList = m_emap[i];
101 const CaloDetDescrElement* calodde = m_cddMgr->get_element(IdentifierHash(i));
102 double eta=calodde->eta();
103 double phi=calodde->phi();
104 for(unsigned int iPart=0;iPart<etaPart.size();++iPart)
105 {
106 double deltaPhi=fmod(phiPart[iPart]-phi+3.0*M_PI,2.0*M_PI)-M_PI;
107 double deltaEta=etaPart[iPart]-eta;
108 if( std::fabs(deltaPhi)<dphi/2. &&
109 std::fabs(deltaEta)<deta/2. )
110 {
111 theLArHitList.setInWindows();
112 break;
113 }
114 } // loop over particles
115 } // loop over cells
116 return true;
117}
#define M_PI
Scalar eta() const
pseudorapidity method
Scalar deltaPhi(const MatrixBase< Derived > &vec) const
Scalar phi() const
phi method
DataModel_detail::const_iterator< DataVector > const_iterator
Definition DataVector.h:838
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
void setInWindows()
Definition LArHitList.h:27
bool is_simulation_particle(const T &p)
Method to establish if a particle (or barcode) was created during the simulation (TODO update to be s...
bool isPhoton(const T &p)
bool isElectron(const T &p)
bool isStable(const T &p)
Identify if the particle is stable, i.e. has not decayed.
double deltaEta(const I4Momentum &p1, const I4Momentum &p2)
Computes efficiently .
Definition P4Helpers.h:66

◆ GetCell()

const LArHitList & LArHitEMap::GetCell ( const unsigned int index) const
inline

Definition at line 43 of file LArHitEMap.h.

43{return m_emap[index];} ;

◆ GetDigit()

const LArDigit * LArHitEMap::GetDigit ( unsigned int index) const
inline

Definition at line 48 of file LArHitEMap.h.

48 {
49 if (index<m_digmap.size()) return m_digmap[index];
50 else return 0;
51 }

◆ GetNbCells()

size_t LArHitEMap::GetNbCells ( void ) const
inline

Definition at line 42 of file LArHitEMap.h.

42{return m_emap.size(); }

◆ GetTimeE()

const std::vector< std::pair< float, float > > & LArHitEMap::GetTimeE ( const IdentifierHash index) const
inline

Definition at line 44 of file LArHitEMap.h.

44{ return m_emap[index].getData();}

Member Data Documentation

◆ m_cabling

const LArOnOffIdMapping* LArHitEMap::m_cabling
private

Definition at line 30 of file LArHitEMap.h.

◆ m_calocell_id

const CaloCell_ID* LArHitEMap::m_calocell_id
private

Definition at line 31 of file LArHitEMap.h.

◆ m_cddMgr

const CaloDetDescrManager* LArHitEMap::m_cddMgr
private

Definition at line 32 of file LArHitEMap.h.

◆ m_digmap

DIGMAP LArHitEMap::m_digmap
private

Definition at line 29 of file LArHitEMap.h.

◆ m_emap

std::vector<LArHitList> LArHitEMap::m_emap
private

Definition at line 28 of file LArHitEMap.h.


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