ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/CmsHi/JetAnalysis/src/HiPFCandAnalyzer.cc
Revision: 1.2
Committed: Tue Sep 20 21:00:45 2011 UTC (13 years, 7 months ago) by yilmaz
Content type: text/plain
Branch: MAIN
Changes since 1.1: +13 -4 lines
Log Message:
PF candidate analyzer

File Contents

# User Rev Content
1 yilmaz 1.1 // -*- C++ -*-
2     //
3     // Package: HiPFCandAnalyzer
4     // Class: HiPFCandAnalyzer
5     //
6     /**\class HiPFCandAnalyzer HiPFCandAnalyzer.cc ana/HiPFCandAnalyzer/src/HiPFCandAnalyzer.cc
7    
8     Description: [one line class summary]
9    
10     Implementation:
11     [Notes on implementation]
12     */
13     //
14     // Original Author: Matt, Nguyen
15     // Created: Oct 10 2010
16     //
17     //
18     //
19    
20    
21     // system include files
22     #include <memory>
23    
24     // stl
25     #include <algorithm>
26    
27     // user include files
28     #include "FWCore/Framework/interface/Frameworkfwd.h"
29     #include "FWCore/Framework/interface/EDAnalyzer.h"
30    
31     #include "FWCore/Framework/interface/Event.h"
32     #include "FWCore/Framework/interface/MakerMacros.h"
33    
34     #include "FWCore/ParameterSet/interface/ParameterSet.h"
35    
36     // ana
37     #include "CmsHi/JetAnalysis/interface/HiPFCandAnalyzer.h"
38    
39     #include "FWCore/ServiceRegistry/interface/Service.h"
40     #include "CommonTools/UtilAlgos/interface/TFileService.h"
41     #include "RecoJets/JetAlgorithms/interface/JetAlgoHelper.h"
42     #include "FWCore/Common/interface/TriggerNames.h"
43     #include "DataFormats/Common/interface/TriggerResults.h"
44     #include "DataFormats/TrackReco/interface/Track.h"
45     #include "DataFormats/VertexReco/interface/Vertex.h"
46     #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
47     #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
48     #include "DataFormats/Math/interface/LorentzVector.h"
49     #include "DataFormats/Math/interface/deltaR.h"
50     #include "DataFormats/Math/interface/deltaPhi.h"
51     #include "TMath.h"
52     #include "TStopwatch.h"
53    
54     #include "DataFormats/HeavyIonEvent/interface/Centrality.h"
55     #include "SimDataFormats/HiGenData/interface/GenHIEvent.h"
56    
57     using namespace std;
58     using namespace edm;
59     using namespace reco;
60    
61     //
62     // constructors and destructor
63     //
64     HiPFCandAnalyzer::HiPFCandAnalyzer(const edm::ParameterSet& iConfig)
65     {
66     // Event source
67     // Event Info
68     pfCandidateLabel_ = iConfig.getParameter<edm::InputTag>("pfCandidateLabel");
69     genLabel_ = iConfig.getParameter<edm::InputTag>("genLabel");
70     jetLabel_ = iConfig.getParameter<edm::InputTag>("jetLabel");
71    
72     pfPtMin_ = iConfig.getParameter<double>("pfPtMin");
73     genPtMin_ = iConfig.getParameter<double>("genPtMin");
74     jetPtMin_ = iConfig.getParameter<double>("jetPtMin");
75    
76     // debug
77     verbosity_ = iConfig.getUntrackedParameter<int>("verbosity", 0);
78    
79 yilmaz 1.2 doJets_ = iConfig.getUntrackedParameter<bool>("doJets",0);
80     doMC_ = iConfig.getUntrackedParameter<bool>("doMC",0);
81 yilmaz 1.1
82     }
83    
84    
85     HiPFCandAnalyzer::~HiPFCandAnalyzer()
86     {
87    
88     // do anything here that needs to be done at desctruction time
89     // (e.g. close files, deallocate resources etc.)
90    
91     }
92    
93    
94     //
95     // member functions
96     //
97    
98     // ------------ method called to for each event ------------
99     void
100     HiPFCandAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
101     {
102    
103     pfEvt_.Clear();
104    
105     // Fill PF info
106    
107     edm::Handle<reco::PFCandidateCollection> pfCandidates;
108     iEvent.getByLabel(pfCandidateLabel_,pfCandidates);
109     const reco::PFCandidateCollection *pfCandidateColl = &(*pfCandidates);
110    
111    
112     for(unsigned icand=0;icand<pfCandidateColl->size(); icand++) {
113     const reco::PFCandidate pfCandidate = pfCandidateColl->at(icand);
114    
115     double pt = pfCandidate.pt();
116    
117     if(pt>pfPtMin_){
118     pfEvt_.pfId_[pfEvt_.nPFpart_] = pfCandidate.particleId();
119     pfEvt_.pfPt_[pfEvt_.nPFpart_] = pt;
120     pfEvt_.pfEta_[pfEvt_.nPFpart_] = pfCandidate.eta();
121     pfEvt_.pfPhi_[pfEvt_.nPFpart_] = pfCandidate.phi();
122     pfEvt_.nPFpart_++;
123     }
124     }
125    
126    
127     // Fill GEN info
128 yilmaz 1.2 if(doMC_){
129 yilmaz 1.1 edm::Handle<reco::GenParticleCollection> genParticles;
130     iEvent.getByLabel(genLabel_,genParticles);
131     const reco::GenParticleCollection* genColl= &(*genParticles);
132    
133     for(unsigned igen=0;igen<genColl->size(); igen++) {
134    
135     const reco::GenParticle gen = genColl->at(igen);
136     double eta = gen.eta();
137     double pt = gen.pt();
138    
139     if(gen.status()==1 && fabs(eta)<3.0 && pt> genPtMin_){
140     pfEvt_.genPDGId_[pfEvt_.nGENpart_] = gen.pdgId();
141     pfEvt_.genPt_[pfEvt_.nGENpart_] = pt;
142     pfEvt_.genEta_[pfEvt_.nGENpart_] = eta;
143     pfEvt_.genPhi_[pfEvt_.nGENpart_] = gen.phi();
144     pfEvt_.nGENpart_++;
145     }
146     }
147 yilmaz 1.2 }
148 yilmaz 1.1
149     // Fill Jet info
150 yilmaz 1.2 if(doJets_){
151 yilmaz 1.1 edm::Handle<pat::JetCollection> jets;
152     iEvent.getByLabel(jetLabel_,jets);
153     const pat::JetCollection *jetColl = &(*jets);
154    
155    
156     for(unsigned ijet=0;ijet<jetColl->size(); ijet++) {
157     const pat::Jet jet = jetColl->at(ijet);
158    
159     double pt = jet.pt();
160    
161     if(pt>jetPtMin_){
162     pfEvt_.pfPt_[pfEvt_.njets_] = pt;
163     pfEvt_.pfEta_[pfEvt_.njets_] = jet.eta();
164     pfEvt_.pfPhi_[pfEvt_.njets_] = jet.phi();
165     pfEvt_.njets_++;
166     }
167     }
168 yilmaz 1.2 }
169 yilmaz 1.1
170     // All done
171     pfTree_->Fill();
172     }
173    
174     /*
175     void HiPFCandAnalyzer::FillEventInfo(const edm::Event& iEvent, const edm::EventSetup& iSetup, TreePFCandEventData & tr)
176     {
177     // General Info
178     tr.run_ = iEvent.id().run();
179     tr.evt_ = iEvent.id().event();
180     tr.lumi_ = iEvent.luminosityBlock();
181    
182     if(!genOnly_&&sampleType_<10){
183     // HI Event info
184     edm::Handle<reco::Centrality> cent;
185     iEvent.getByLabel(edm::InputTag("hiCentrality"),cent);
186     Double_t hf = cent->EtHFhitSum();
187     // Get Centrality bin
188     cbins_ = getCentralityBinsFromDB(iSetup);
189     tr.cent_ = cbins_->getBin(hf)*(100./cbins_->getNbins());
190     }
191    
192     if (isMC_&&sampleType_<10) {
193     edm::Handle<edm::GenHIEvent> mchievt;
194     iEvent.getByLabel(edm::InputTag("heavyIon"),mchievt);
195     tr.b_ = mchievt->b();
196     tr.npart_ = mchievt->Npart();
197     tr.ncoll_ = mchievt->Ncoll();
198     }
199     }
200     */
201     void HiPFCandAnalyzer::beginJob()
202     {
203    
204     // -- trees --
205     pfTree_ = fs->make<TTree>("pfTree","dijet tree");
206     pfEvt_.SetTree(pfTree_);
207 yilmaz 1.2 pfEvt_.doMC = doMC_;
208     pfEvt_.doJets = doJets_;
209    
210 yilmaz 1.1 pfEvt_.SetBranches();
211     }
212    
213     // ------------ method called once each job just after ending the event loop ------------
214     void
215     HiPFCandAnalyzer::endJob() {
216     // ===== Done =====
217     /* if (verbosity_>=1) {
218     cout << endl << "================ Ana Process Summaries =============" << endl;
219     cout << " AnaJet: " << jetsrc_ << endl;
220     if (refJetType_>=0) cout << " RefJet: " << refjetsrc_ << endl;
221     cout << " AnaTrk: " << trksrc_ << endl;
222     cout << "# HI Events : "<< numHiEvtSel_<< endl;
223     cout << "# Base Events: "<< numEvtSel_ << endl;
224     cout << "# Jet Events: "<< numJetEvtSel_<< endl;
225     }
226     */
227     }
228    
229     // constructors
230     TreePFCandEventData::TreePFCandEventData(){
231    
232     }
233    
234    
235     // set branches
236     void TreePFCandEventData::SetBranches()
237     {
238     // --event level--
239    
240     // -- particle info --
241     tree_->Branch("nPFpart",&(this->nPFpart_),"nPFpart/I");
242     tree_->Branch("pfId",this->pfId_,"pfId[nPFpart]/I");
243     tree_->Branch("pfPt",this->pfPt_,"pfPt[nPFpart]/F");
244     tree_->Branch("pfEta",this->pfEta_,"pfEta[nPFpart]/F");
245     tree_->Branch("pfPhi",this->pfPhi_,"pfPhi[nPFpart]/F");
246    
247     // -- jet info --
248 yilmaz 1.2 if(doJets){
249 yilmaz 1.1 tree_->Branch("njets",&(this->njets_),"njets/I");
250     tree_->Branch("jetPt",this->jetPt_,"jetPt[njets]/F");
251     tree_->Branch("jetEta",this->jetEta_,"jetEta[njets]/F");
252     tree_->Branch("jetPhi",this->jetPhi_,"jetPhi[njets]/F");
253 yilmaz 1.2 }
254 yilmaz 1.1
255     // -- gen info --
256 yilmaz 1.2 if(doMC){
257 yilmaz 1.1 tree_->Branch("nGENpart",&(this->nGENpart_),"nGENpart/I");
258     tree_->Branch("genPDGId",this->genPDGId_,"genPDGId[nGENpart]/I");
259     tree_->Branch("genPt",this->genPt_,"genPt[nGENpart]/F");
260     tree_->Branch("genEta",this->genEta_,"genEta[nGENpart]/F");
261     tree_->Branch("genPhi",this->genPhi_,"genPhi[nGENpart]/F");
262 yilmaz 1.2 }
263 yilmaz 1.1
264     }
265     void TreePFCandEventData::Clear()
266     {
267     // event
268    
269     nPFpart_ = 0;
270     njets_ = 0;
271     nGENpart_ = 0;
272     }
273    
274     DEFINE_FWK_MODULE(HiPFCandAnalyzer);