CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
LhcTrackAnalyzer.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: LhcTrackAnalyzer
4 // Class: LhcTrackAnalyzer
5 //
16 //
17 
18 // updated to 25/2/2009 5.30 pm
19 
20 //
21 //
22 
23 // system include files
24 #include <memory>
25 
26 
27 // user include files
31 
36 
38 
39 
40 #include "TH1F.h"
41 #include "TH2F.h"
42 #include "TFile.h"
43 #include "TROOT.h"
44 #include "TChain.h"
45 #include "TNtuple.h"
55 
60 
61 // Constructor
62 
64 
65 {
66  //now do what ever initialization is needed
67  debug_ = iConfig.getParameter<bool> ("Debug");
68  TrackCollectionTag_ = iConfig.getParameter<edm::InputTag>("TrackCollectionTag");
69  PVtxCollectionTag_ = iConfig.getParameter<edm::InputTag>("PVtxCollectionTag");
70  filename_ = iConfig.getParameter<std::string>("OutputFileName");
71 }
72 
73 // Destructor
75 {
76 
77  // do anything here that needs to be done at desctruction time
78  // (e.g. close files, deallocate resources etc.)
79 
80 }
81 
82 
83 //
84 // member functions
85 //
86 
87 // ------------ method called to for each event ------------
88 void
90 {
91  using namespace edm;
92  using namespace reco;
93  using namespace std;
94 
95  //=======================================================
96  // Initialize Root-tuple variables
97  //=======================================================
98 
99  SetVarToZero();
100 
101  //=======================================================
102  // Retrieve the Track information
103  //=======================================================
104 
105  Handle< TrackCollection> trackCollectionHandle;
106  iEvent.getByLabel(TrackCollectionTag_, trackCollectionHandle);
107  Handle<VertexCollection> vertexCollectionHandle;
108  iEvent.getByLabel(PVtxCollectionTag_, vertexCollectionHandle);
109  for(VertexCollection::const_iterator vtx = vertexCollectionHandle->begin();vtx!=vertexCollectionHandle->end(); ++vtx)
110  {
111  if(vtx==vertexCollectionHandle->begin()){
112  if(vtx->isFake())goodvtx_=false;
113  else goodvtx_=true;
114  }
115  else break;
116  }
117 
118 
119 
120  goodbx_=true;
121  // int bx = iEvent.bunchCrossing();
122  //if (bx==51 || bx==2724) goodbx_=true;
123 
124 
125 
126  run_=iEvent.id().run();
127  event_=iEvent.id().event();
128 
129  if(debug_)
130  cout<<"LhcTrackAnalyzer::analyze() looping over "<< trackCollectionHandle->size()<< "tracks." << endl;
131 
132  // unsigned int i = 0;
133  bool toomanytracks=false;
134  for(TrackCollection::const_iterator track = trackCollectionHandle->begin(); track!= trackCollectionHandle->end(); ++track)
135  {
136  //if(!toomanytracks){
137  if ( nTracks_ >= nMaxtracks_ ) {
138  std::cout << " LhcTrackAnalyzer::analyze() : Warning - Run "<< run_<<" Event "<< event_<<"\tNumber of tracks: " << trackCollectionHandle->size()<< " , greater than " << nMaxtracks_ << std::endl;
139  toomanytracks=true;
140  continue;
141 
142  }
143  // else{
144  pt_[nTracks_] = track->pt();
145  eta_[nTracks_] = track->eta();
146  phi_[nTracks_] = track->phi();
147  chi2_[nTracks_] = track->chi2();
148  chi2ndof_[nTracks_] = track->normalizedChi2();
149  charge_[nTracks_] = track->charge();
150  qoverp_[nTracks_] = track->qoverp();
151  dz_[nTracks_] = track->dz();
152  dxy_[nTracks_] = track->dxy();
153  xPCA_[nTracks_] = track->vertex().x();
154  yPCA_[nTracks_] = track->vertex().y();
155  zPCA_[nTracks_] = track->vertex().z();
156  validhits_[nTracks_][0]=track->numberOfValidHits();
157  validhits_[nTracks_][1]=track->hitPattern().numberOfValidPixelBarrelHits();
158  validhits_[nTracks_][2]=track->hitPattern().numberOfValidPixelEndcapHits();
159  validhits_[nTracks_][3]=track->hitPattern().numberOfValidStripTIBHits();
160  validhits_[nTracks_][4]=track->hitPattern().numberOfValidStripTIDHits();
161  validhits_[nTracks_][5]=track->hitPattern().numberOfValidStripTOBHits();
162  validhits_[nTracks_][6]=track->hitPattern().numberOfValidStripTECHits();
163 
164 
165 
166  int myalgo=-88;
167  if(track->algo()==reco::TrackBase::undefAlgorithm)myalgo=0;
168  if(track->algo()==reco::TrackBase::ctf)myalgo=1;
169  if(track->algo()==reco::TrackBase::iter0)myalgo=4;
170  if(track->algo()==reco::TrackBase::iter1)myalgo=5;
171  if(track->algo()==reco::TrackBase::iter2)myalgo=6;
172  if(track->algo()==reco::TrackBase::iter3)myalgo=7;
173  if(track->algo()==reco::TrackBase::iter4)myalgo=8;
174  if(track->algo()==reco::TrackBase::iter5)myalgo=9;
175  if(track->algo()==reco::TrackBase::iter6)myalgo=10;
176  if(track->algo()==reco::TrackBase::iter7)myalgo=11;
177  trkAlgo_[nTracks_] = myalgo;
178 
179  int myquality=-99;
180  if(track->quality(reco::TrackBase::undefQuality))myquality=-1;
181  if(track->quality(reco::TrackBase::loose))myquality=0;
182  if(track->quality(reco::TrackBase::tight))myquality=1;
183  if(track->quality(reco::TrackBase::highPurity))myquality=2;
184  //if(track->quality(reco::TrackBase::confirmed))myquality=3;
185  // if(track->quality(reco::TrackBase::goodIterative))myquality=4;
186  // if(track->quality(reco::TrackBase::qualitySize))myquality=5;
187  trkQuality_[nTracks_]= myquality;
188 
189  if(track->quality(reco::TrackBase::highPurity))isHighPurity_[nTracks_]=1;
190  else isHighPurity_[nTracks_]=0;
191  nTracks_++;
192 
193 
194 
195 
196  // }//end if not too many tracks
197 
198 
199 
200  }//end loop on tracks
201 
202  for(int d=0;d<nTracks_;++d){
203  if(abs(trkQuality_[d])>5)cout<<"MYQUALITY!!! " <<trkQuality_[d] <<" at track # "<<d<<"/"<< nTracks_<<endl;
204  }
205 
206 
207 
208 
209  rootTree_->Fill();
210 }
211 
212 
213 // ------------ method called once each job before begining the event loop ------------
215 {
216  edm::LogInfo("beginJob") << "Begin Job" << std::endl;
217  // Define TTree for output
218  rootFile_ = new TFile(filename_.c_str(),"recreate");
219  rootTree_ = new TTree("tree","Lhc Track tree");
220 
221  // Track Paramters
222  rootTree_->Branch("run",&run_,"run/I");
223  rootTree_->Branch("event",&event_,"event/I");
224  rootTree_->Branch("goodbx",&goodbx_,"goodbx/O");
225  rootTree_->Branch("goodvtx",&goodvtx_,"goodvtx/O");
226  rootTree_->Branch("nTracks",&nTracks_,"nTracks/I");
227  rootTree_->Branch("pt",&pt_,"pt[nTracks]/D");
228  rootTree_->Branch("eta",&eta_,"eta[nTracks]/D");
229  rootTree_->Branch("phi",&phi_,"phi[nTracks]/D");
230  rootTree_->Branch("chi2",&chi2_,"chi2[nTracks]/D");
231  rootTree_->Branch("chi2ndof",&chi2ndof_,"chi2ndof[nTracks]/D");
232  rootTree_->Branch("charge",&charge_,"charge[nTracks]/I");
233  rootTree_->Branch("qoverp",&qoverp_,"qoverp[nTracks]/D");
234  rootTree_->Branch("dz",&dz_,"dz[nTracks]/D");
235  rootTree_->Branch("dxy",&dxy_,"dxy[nTracks]/D");
236  rootTree_->Branch("xPCA",&xPCA_,"xPCA[nTracks]/D");
237  rootTree_->Branch("yPCA",&yPCA_,"yPCA[nTracks]/D");
238  rootTree_->Branch("zPCA",&zPCA_,"zPCA[nTracks]/D");
239  rootTree_->Branch("isHighPurity",&isHighPurity_,"isHighPurity[nTracks]/I");
240  rootTree_->Branch("trkQuality",&trkQuality_,"trkQuality[nTracks]/I");
241  rootTree_->Branch("trkAlgo",&trkAlgo_,"trkAlgo[nTracks]/I");
242  rootTree_->Branch("nValidHits",&validhits_,"nValidHits[nTracks][7]/I");
243 
244 
245 }
246 
247 // ------------ method called once each job just after ending the event loop ------------
249 {
250  if ( rootFile_ ) {
251  rootFile_->Write();
252  rootFile_->Close();
253  }
254 
255 
256 
257 }
258 
260  run_=-1;
261  event_=-99;
262  nTracks_ = 0;
263  for ( int i=0; i<nMaxtracks_; ++i ) {
264  pt_[i] = 0;
265  eta_[i] = 0;
266  phi_[i] = 0;
267  chi2_[i] = 0;
268  chi2ndof_[i] = 0;
269  charge_[i] = 0;
270  qoverp_[i] = 0;
271  dz_[i] = 0;
272  dxy_[i] = 0;
273  xPCA_[i] = 0;
274  yPCA_[i] = 0;
275  zPCA_[i] = 0;
276  trkQuality_[i] = 0;
277  trkAlgo_[i] = -1;
278  isHighPurity_[i]=-3;
279  for(int j=0;j<7;j++){
280  validhits_[nTracks_][j]=-1*j;
281  }
282  }
283 
284 
285 }
286 
287 //define this as a plug-in
RunNumber_t run() const
Definition: EventID.h:42
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:44
int i
Definition: DBlmapReader.cc:9
int validhits_[nMaxtracks_][7]
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
#define abs(x)
Definition: mlp_lapack.h:159
int isHighPurity_[nMaxtracks_]
double dz_[nMaxtracks_]
double zPCA_[nMaxtracks_]
edm::InputTag TrackCollectionTag_
int trkQuality_[nMaxtracks_]
double qoverp_[nMaxtracks_]
virtual void endJob()
double eta_[nMaxtracks_]
int iEvent
Definition: GenABIO.cc:243
virtual void beginJob()
std::string filename_
double phi_[nMaxtracks_]
double chi2ndof_[nMaxtracks_]
double yPCA_[nMaxtracks_]
int j
Definition: DBlmapReader.cc:9
double dxy_[nMaxtracks_]
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:355
int trkAlgo_[nMaxtracks_]
edm::InputTag PVtxCollectionTag_
virtual void analyze(const edm::Event &, const edm::EventSetup &)
double xPCA_[nMaxtracks_]
edm::EventID id() const
Definition: EventBase.h:56
static const int nMaxtracks_
double chi2_[nMaxtracks_]
int charge_[nMaxtracks_]
tuple cout
Definition: gather_cfg.py:41
LhcTrackAnalyzer(const edm::ParameterSet &)
double pt_[nMaxtracks_]