CMS 3D CMS Logo
List of all members | Public Member Functions | Private Member Functions | Private Attributes
ElectronMCTruthFinder Class Reference

#include <ElectronMCTruthFinder.h>

Public Member Functions

 ElectronMCTruthFinder ()
 
std::vector< ElectronMCTruth > find (const std::vector< SimTrack > &simTracks, const std::vector< SimVertex > &simVertices)
 
virtual ~ElectronMCTruthFinder ()
 

Private Member Functions

void fill (const std::vector< SimTrack > &theSimTracks, const std::vector< SimVertex > &theSimVertices)
 

Private Attributes

std::map< unsigned, unsigned > geantToIndex_
 
PhotonMCTruthFinder * thePhotonMCTruthFinder_
 

Detailed Description

Definition at line 13 of file ElectronMCTruthFinder.h.

Constructor & Destructor Documentation

◆ ElectronMCTruthFinder()

ElectronMCTruthFinder::ElectronMCTruthFinder ( )

Definition at line 5 of file ElectronMCTruthFinder.cc.

5 {}

◆ ~ElectronMCTruthFinder()

ElectronMCTruthFinder::~ElectronMCTruthFinder ( )
virtual

Definition at line 7 of file ElectronMCTruthFinder.cc.

7 {}

Member Function Documentation

◆ fill()

void ElectronMCTruthFinder::fill ( const std::vector< SimTrack > &  theSimTracks,
const std::vector< SimVertex > &  theSimVertices 
)
private

Definition at line 195 of file ElectronMCTruthFinder.cc.

References geantToIndex_, L1TkPrimaryVertex_cfi::nVtx, TrackCandidateProducer_cfi::simTracks, and HGCalValidator_cfi::simVertices.

Referenced by find().

195  {
196  //std::cout << " ElectronMCTruthFinder::fill " << std::endl;
197 
198  unsigned nVtx = simVertices.size();
199  unsigned nTks = simTracks.size();
200 
201  // Empty event, do nothin'
202  if (nVtx == 0)
203  return;
204 
205  // create a map associating geant particle id and position in the
206  // event SimTrack vector
207  for (unsigned it = 0; it < nTks; ++it) {
208  geantToIndex_[simTracks[it].trackId()] = it;
209  //std::cout << " ElectronMCTruthFinder::fill it " << it << " simTracks[it].trackId() " << simTracks[it].trackId() << std::endl;
210  }
211 }
ElectronMCTruthFinder::geantToIndex_
std::map< unsigned, unsigned > geantToIndex_
Definition: ElectronMCTruthFinder.h:24
L1TkPrimaryVertex_cfi.nVtx
nVtx
Definition: L1TkPrimaryVertex_cfi.py:19
TrackCandidateProducer_cfi.simTracks
simTracks
Definition: TrackCandidateProducer_cfi.py:15
HGCalValidator_cfi.simVertices
simVertices
Definition: HGCalValidator_cfi.py:66

◆ find()

std::vector< ElectronMCTruth > ElectronMCTruthFinder::find ( const std::vector< SimTrack > &  simTracks,
const std::vector< SimVertex > &  simVertices 
)

Now store the electron truth

here fill the electron

Definition at line 9 of file ElectronMCTruthFinder.cc.

References funct::abs(), MillePedeFileConverter_cfg::e, simKBmtfDigis_cfi::eLoss, fill(), geantToIndex_, if(), CoreSimTrack::momentum(), objects.autophobj::motherId, SimVertex::parentIndex(), CoreSimVertex::position(), mps_fire::result, CoreSimTrack::trackId(), CoreSimTrack::type(), and SimTrack::vertIndex().

Referenced by MCElectronAnalyzer::analyze().

10  {
11  //std::cout << " ElectronMCTruthFinder::find " << std::endl;
12 
13  std::vector<ElectronMCTruth> result;
14 
15  // Local variables
16  //const int SINGLE=1;
17  //const int DOUBLE=2;
18  //const int PYTHIA=3;
19  //const int ELECTRON_FLAV=1;
20  //const int PIZERO_FLAV=2;
21  //const int PHOTON_FLAV=3;
22 
23  //int ievtype=0;
24  //int ievflav=0;
25 
26  std::vector<SimTrack> electronTracks;
27  SimVertex primVtx;
28 
29  fill(theSimTracks, theSimVertices);
30 
31  int iPV = -1;
32  //int partType1=0;
33  //int partType2=0;
34  std::vector<SimTrack>::const_iterator iFirstSimTk = theSimTracks.begin();
35  if (!(*iFirstSimTk).noVertex()) {
36  iPV = (*iFirstSimTk).vertIndex();
37 
38  int vtxId = (*iFirstSimTk).vertIndex();
39  primVtx = theSimVertices[vtxId];
40 
41  //partType1 = (*iFirstSimTk).type();
42 
43  //std::cout << " Primary vertex id " << iPV << " first track type " << (*iFirstSimTk).type() << std::endl;
44  } else {
45  //std::cout << " First track has no vertex " << std::endl;
46  }
47 
48  // CLHEP::HepLorentzVector primVtxPos= primVtx.position();
49  math::XYZTLorentzVectorD primVtxPos(
50  primVtx.position().x(), primVtx.position().y(), primVtx.position().z(), primVtx.position().e());
51 
52  // Look at a second track
53  iFirstSimTk++;
54  //if ( iFirstSimTk!= theSimTracks.end() ) {
55  //
56  // if ( (*iFirstSimTk).vertIndex() == iPV) {
57  // partType2 = (*iFirstSimTk).type();
58  //
59  // //std::cout << " second track type " << (*iFirstSimTk).type() << " vertex " << (*iFirstSimTk).vertIndex() << std::endl;
60  //
61  // } else {
62  // //std::cout << " Only one kine track from Primary Vertex " << std::endl;
63  // }
64  //}
65 
66  //std::cout << " Loop over all particles " << std::endl;
67 
68  int npv = 0;
69  for (std::vector<SimTrack>::const_iterator iSimTk = theSimTracks.begin(); iSimTk != theSimTracks.end(); ++iSimTk) {
70  if ((*iSimTk).noVertex())
71  continue;
72 
73  //int vertexId = (*iSimTk).vertIndex();
74  //SimVertex vertex = theSimVertices[vertexId];
75 
76  //std::cout << " Particle type " << (*iSimTk).type() << " Sim Track ID " << (*iSimTk).trackId() << " momentum " << (*iSimTk).momentum() << " vertex position " << vertex.position() << " vertex ID " << vertexId << std::endl;
77  if ((*iSimTk).vertIndex() == iPV) {
78  npv++;
79  if (std::abs((*iSimTk).type()) == 11) {
80  //std::cout << " Found a primary electron with ID " << (*iSimTk).trackId() << " momentum " << (*iSimTk).momentum() << std::endl;
81  electronTracks.push_back(*iSimTk);
82  }
83  }
84  }
85  //std::cout << " There are " << npv << " particles originating in the PV " << std::endl;
86 
87  //if(npv > 4) {
88  // ievtype = PYTHIA;
89  //} else if(npv == 1) {
90  // if( std::abs(partType1) == 11 ) {
91  // ievtype = SINGLE; ievflav = ELECTRON_FLAV;
92  // } else if(partType1 == 111) {
93  // ievtype = SINGLE; ievflav = PIZERO_FLAV;
94  // } else if(partType1 == 22) {
95  // ievtype = SINGLE; ievflav = PHOTON_FLAV;
96  // }
97  //} else if(npv == 2) {
98  // if ( std::abs(partType1) == 11 && std::abs(partType2) == 11 ) {
99  // ievtype = DOUBLE; ievflav = ELECTRON_FLAV;
100  // } else if(partType1 == 111 && partType2 == 111) {
101  // ievtype = DOUBLE; ievflav = PIZERO_FLAV;
102  // } else if(partType1 == 22 && partType2 == 22) {
103  // ievtype = DOUBLE; ievflav = PHOTON_FLAV;
104  // }
105  //}
106 
108  std::vector<CLHEP::Hep3Vector> bremPos;
109  std::vector<CLHEP::HepLorentzVector> pBrem;
110  std::vector<float> xBrem;
111 
112  for (std::vector<SimTrack>::iterator iEleTk = electronTracks.begin(); iEleTk != electronTracks.end(); ++iEleTk) {
113  //std::cout << " Looping on the primary electron pt " << std::sqrt((*iEleTk).momentum().perp2()) << " electron track ID " << (*iEleTk).trackId() << std::endl;
114 
115  SimTrack trLast = (*iEleTk);
116  unsigned int eleId = (*iEleTk).trackId();
117  float remainingEnergy = trLast.momentum().e();
118  // CLHEP::HepLorentzVector motherMomentum = (*iEleTk).momentum();
119  // CLHEP::HepLorentzVector primEleMom = (*iEleTk).momentum();
120  math::XYZTLorentzVectorD motherMomentum(
121  (*iEleTk).momentum().x(), (*iEleTk).momentum().y(), (*iEleTk).momentum().z(), (*iEleTk).momentum().e());
122  math::XYZTLorentzVectorD primEleMom(motherMomentum);
123  int eleVtxIndex = (*iEleTk).vertIndex();
124 
125  bremPos.clear();
126  pBrem.clear();
127  xBrem.clear();
128 
129  for (std::vector<SimTrack>::const_iterator iSimTk = theSimTracks.begin(); iSimTk != theSimTracks.end(); ++iSimTk) {
130  if ((*iSimTk).noVertex())
131  continue;
132  if ((*iSimTk).vertIndex() == iPV)
133  continue;
134 
135  //std::cout << " (*iEleTk)->trackId() " << (*iEleTk).trackId() << " (*iEleTk)->vertIndex() "<< (*iEleTk).vertIndex() << " (*iSimTk).vertIndex() " << (*iSimTk).vertIndex() << " (*iSimTk).type() " << (*iSimTk).type() << " (*iSimTk).trackId() " << (*iSimTk).trackId() << std::endl;
136 
137  int vertexId1 = (*iSimTk).vertIndex();
138  SimVertex vertex1 = theSimVertices[vertexId1];
139  int vertexId2 = trLast.vertIndex();
140  //SimVertex vertex2 = theSimVertices[vertexId2];
141 
142  int motherId = -1;
143 
144  if ((vertexId1 == vertexId2) && ((*iSimTk).type() == (*iEleTk).type()) && trLast.type() == 22) {
145  //std::cout << " Here a e/gamma brem vertex " << std::endl;
146 
147  //std::cout << " Secondary from electron: particle1 type " << (*iSimTk).type() << " trackId " << (*iSimTk).trackId() << " vertex ID " << vertexId1 << " vertex position " << std::sqrt(vertex1.position().perp2()) << " parent index "<< vertex1.parentIndex() << std::endl;
148 
149  //std::cout << " Secondary from electron: particle2 type " << trLast.type() << " trackId " << trLast.trackId()<< " vertex ID " << vertexId2 << " vertex position " << std::sqrt(vertex2.position().perp2()) << " parent index " << vertex2.parentIndex() << std::endl;
150 
151  //std::cout << " Electron pt " << std::sqrt((*iSimTk).momentum().perp2()) << " photon pt " << std::sqrt(trLast.momentum().perp2()) << "Mother electron pt " << sqrt(motherMomentum.perp2()) << std::endl;
152  //std::cout << " eleId " << eleId << std::endl;
153  float eLoss = remainingEnergy - ((*iSimTk).momentum() + trLast.momentum()).e();
154  //std::cout << " eLoss " << eLoss << std::endl;
155 
156  if (vertex1.parentIndex()) {
157  unsigned motherGeantId = vertex1.parentIndex();
158  std::map<unsigned, unsigned>::iterator association = geantToIndex_.find(motherGeantId);
159  if (association != geantToIndex_.end())
160  motherId = association->second;
161 
162  //int motherType = motherId == -1 ? 0 : theSimTracks[motherId].type();
163  //std::cout << " Parent to this vertex motherId " << motherId << " mother type " << motherType << " Sim track ID " << theSimTracks[motherId].trackId() << std::endl;
164  if (theSimTracks[motherId].trackId() == eleId) {
165  //std::cout << " ***** Found the Initial Mother Electron **** theSimTracks[motherId].trackId() " << theSimTracks[motherId].trackId() << " eleId " << eleId << std::endl;
166  eleId = (*iSimTk).trackId();
167  remainingEnergy = (*iSimTk).momentum().e();
168  motherMomentum = (*iSimTk).momentum();
169 
170  pBrem.push_back(CLHEP::HepLorentzVector(
171  trLast.momentum().px(), trLast.momentum().py(), trLast.momentum().pz(), trLast.momentum().e()));
172  bremPos.push_back(CLHEP::HepLorentzVector(
173  vertex1.position().x(), vertex1.position().y(), vertex1.position().z(), vertex1.position().t()));
174  xBrem.push_back(eLoss);
175  }
176 
177  } else {
178  //std::cout << " This vertex has no parent tracks " << std::endl;
179  }
180  }
181  trLast = (*iSimTk);
182 
183  } // End loop over all SimTracks
184  //std::cout << " Going to build the ElectronMCTruth: pBrem size " << pBrem.size() << std::endl;
186  CLHEP::HepLorentzVector tmpEleMom(primEleMom.px(), primEleMom.py(), primEleMom.pz(), primEleMom.e());
187  CLHEP::HepLorentzVector tmpVtxPos(primVtxPos.x(), primVtxPos.y(), primVtxPos.z(), primVtxPos.t());
188  result.push_back(ElectronMCTruth(tmpEleMom, eleVtxIndex, bremPos, pBrem, xBrem, tmpVtxPos, (*iEleTk)));
189 
190  } // End loop over primary electrons
191 
192  return result;
193 }
ElectronMCTruthFinder::geantToIndex_
std::map< unsigned, unsigned > geantToIndex_
Definition: ElectronMCTruthFinder.h:24
math::XYZTLorentzVectorD
ROOT::Math::LorentzVector< ROOT::Math::PxPyPzE4D< double > > XYZTLorentzVectorD
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:14
mps_fire.result
result
Definition: mps_fire.py:311
CoreSimVertex::position
const math::XYZTLorentzVectorD & position() const
Definition: CoreSimVertex.h:21
ElectronMCTruthFinder::fill
void fill(const std::vector< SimTrack > &theSimTracks, const std::vector< SimVertex > &theSimVertices)
Definition: ElectronMCTruthFinder.cc:195
simKBmtfDigis_cfi.eLoss
eLoss
Definition: simKBmtfDigis_cfi.py:25
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
CoreSimTrack::momentum
const math::XYZTLorentzVectorD & momentum() const
Definition: CoreSimTrack.h:19
CoreSimTrack::type
int type() const
particle type (HEP PDT convension)
Definition: CoreSimTrack.h:22
SimTrack::vertIndex
int vertIndex() const
index of the vertex in the Event container (-1 if no vertex)
Definition: SimTrack.h:33
hgcal::association
std::tuple< layerClusterToCaloParticle, caloParticleToLayerCluster > association
Definition: LCToCPAssociatorByEnergyScoreImpl.h:62
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
CoreSimTrack::trackId
unsigned int trackId() const
Definition: CoreSimTrack.h:31
SimVertex::parentIndex
int parentIndex() const
Definition: SimVertex.h:29

Member Data Documentation

◆ geantToIndex_

std::map<unsigned, unsigned> ElectronMCTruthFinder::geantToIndex_
private

Definition at line 24 of file ElectronMCTruthFinder.h.

Referenced by fill(), and find().

◆ thePhotonMCTruthFinder_

PhotonMCTruthFinder* ElectronMCTruthFinder::thePhotonMCTruthFinder_
private

Definition at line 25 of file ElectronMCTruthFinder.h.

Generated for CMSSW Reference Manual by doxygen 1.8.14