ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MPIAnalyzer/src/MPIntuple.cc
(Generate patch)

Comparing UserCode/MPIAnalyzer/src/MPIntuple.cc (file contents):
Revision 1.1 by naodell, Fri May 21 18:20:49 2010 UTC vs.
Revision 1.4 by naodell, Sat May 22 23:34:09 2010 UTC

# Line 7 | Line 7
7   //
8   // Original Author:  "Nathaniel Odell"
9   // Secondary Author: Steven Won
10 + // With contributions from: Andrey Pozdnyakov
11   //         Created:  Thurs April 22  2010
12   // $Id$
13   //
# Line 30 | Line 31
31   #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
32   #include "Geometry/Records/interface/CaloGeometryRecord.h"
33  
34 + // Jet and vertex functions
35   #include "DataFormats/JetReco/interface/CaloJet.h"
36   #include "DataFormats/JetReco/interface/CaloJetCollection.h"
37   #include "DataFormats/JetReco/interface/PFJet.h"
# Line 40 | Line 42
42   #include "DataFormats/VertexReco/interface/Vertex.h"
43   #include "DataFormats/VertexReco/interface/VertexFwd.h"
44  
45 + // Need to get correctors
46 + #include "JetMETCorrections/Objects/interface/JetCorrector.h"
47 +
48 + // ntuple storage classes
49 + #include "TCJet.h"
50 + #include "TCPrimaryVtx.h"
51 +
52 + // Need for HLT trigger info:
53 + #include "FWCore/Common/interface/TriggerNames.h"
54 + #include "DataFormats/Common/interface/TriggerResults.h"
55 + #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
56 +
57 + //Root  stuff
58   #include "TROOT.h"
59   #include "TFile.h"
60   #include "TTree.h"
# Line 70 | Line 85 | class MPIntuple : public edm::EDAnalyzer
85        virtual void analyze(const edm::Event&, const edm::EventSetup&);
86        virtual void endJob() ;
87  
88 +  bool triggerDecision(edm::Handle<edm::TriggerResults> &hltR, int iTrigger);
89 +
90        // ----------member data ---------------------------
91  
92    edm::InputTag PFJetHandle_;
76  edm::InputTag CaloJetHandle_;
77  edm::InputTag PFJetHandle_corr;
78  edm::InputTag CaloJetHandle_corr;
93    edm::InputTag GenJetHandle_;
94    edm::InputTag PrimaryVtxHandle_;
95 +  edm::InputTag triggerResultsTag_;
96  
97  
98    //  Counting variables   //
# Line 87 | Line 102 | class MPIntuple : public edm::EDAnalyzer
102    TTree *sTree;
103    TFile* ntupleFile;
104  
105 <  TClonesArray* jetP4_ak5PF;
91 <  TClonesArray* jetP4_ak5PF_corr;
92 <  TClonesArray* jetVtx3_ak5PF;
93 <  TClonesArray* jetVtx3_ak5PF_corr;
94 <  TClonesArray* jetP4_ak5Calo;
95 <  TClonesArray* jetP4_ak5Calo_corr;
96 <  TClonesArray* jetVtx3_ak5Calo;
97 <  TClonesArray* jetVtx3_ak5Calo_corr;
105 >  TClonesArray* jet_ak5PF;
106    TClonesArray* jetP4_ak5Gen;
107 <  TClonesArray* jetVtx3_ak5Gen;
100 <  TClonesArray* PrimaryVtx3;
107 >  TClonesArray* vtx_offline;
108  
109    bool doGenJets_;
103  bool doCaloJets_;
110    bool doPFJets_;
111  
112    string rootfilename;
113  
114 +  //Triggers
115 +  string hlTriggerResults_, hltName_, triggerName_;
116 +  TriggerNames triggerNames;
117 +  HLTConfigProvider hltConfig_;
118 +  vector<string>  hlNames;
119 +
120 +
121   };
122  
123   MPIntuple::MPIntuple(const edm::ParameterSet& iConfig):
124  
125    PFJetHandle_(iConfig.getUntrackedParameter<edm::InputTag>("PFJetTag")),
113  CaloJetHandle_(iConfig.getUntrackedParameter<edm::InputTag>("CaloJetTag")),
114  PFJetHandle_corr(iConfig.getUntrackedParameter<edm::InputTag>("PFJetTag_corr")),
115  CaloJetHandle_corr(iConfig.getUntrackedParameter<edm::InputTag>("CaloJetTag_corr")),
126    GenJetHandle_(iConfig.getUntrackedParameter<edm::InputTag>("GenJetTag")),
127    PrimaryVtxHandle_(iConfig.getUntrackedParameter<edm::InputTag>("PrimaryVtxTag")),
128    doGenJets_(iConfig.getUntrackedParameter<bool>("doGenJets")),
119  doCaloJets_(iConfig.getUntrackedParameter<bool>("doCaloJets")),
129    doPFJets_(iConfig.getUntrackedParameter<bool>("doPFJets")),
130 <  rootfilename(iConfig.getUntrackedParameter<string>("rootfilename"))
131 <
130 >  rootfilename(iConfig.getUntrackedParameter<string>("rootfilename")),
131 >  hlTriggerResults_(iConfig.getUntrackedParameter<string>("HLTriggerResults","TriggerResults")),
132 >  hltName_(iConfig.getUntrackedParameter<string>("hltName","HLT"))
133   {
134 +  //edm::TriggerNames
135 +  // triggerNames(iConfig);
136  
137   }
138  
# Line 130 | Line 142 | MPIntuple::~MPIntuple()
142  
143   }
144  
133
145   //
146   // member functions
147   //
# Line 144 | Line 155 | void MPIntuple::analyze(const edm::Event
155    lumiSection  = (unsigned int)iEvent.getLuminosityBlock().luminosityBlock();
156    
157    int PFcount = 0;
147  int Calocount = 0;
148  int PFcount_corr = 0;
149  int Calocount_corr = 0;
158    int Gencount = 0;
159    int Vtxcount = 0;
160  
161    if(doPFJets_){
162 <    
162 >
163 >    const JetCorrector* correctorL2 = JetCorrector::getJetCorrector ("ak5PFL2Relative",iSetup);
164 >    const JetCorrector* correctorL3 = JetCorrector::getJetCorrector ("ak5PFL3Absolute",iSetup);
165 >
166      Handle<reco::PFJetCollection> PFJets;
167      iEvent.getByLabel(PFJetHandle_, PFJets);
168      
# Line 160 | Line 171 | void MPIntuple::analyze(const edm::Event
171      for(PFJetCollection::const_iterator jet_iter = PFJets->begin(); jet_iter!= PFJets->end(); ++jet_iter){
172        
173        reco::PFJet myJet = reco::PFJet(*jet_iter);
174 <      
175 <      if(myJet.et() > 5){
176 <        
177 <        new ((*jetP4_ak5PF)[PFcount]) TLorentzVector(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
167 <        new ((*jetVtx3_ak5PF)[PFcount]) TVector3(myJet.vx(), myJet.vy(), myJet.vz());
168 <        
169 <        ++PFcount;
170 <      }
171 <      
172 <    }
173 <    
174 <    
175 <    PFcount_corr = 0;
176 <    
177 <    Handle<reco::PFJetCollection> PFJets_corr;
178 <    iEvent.getByLabel(PFJetHandle_corr, PFJets_corr);
179 <    
180 <    for(PFJetCollection::const_iterator jet_iter = PFJets_corr->begin(); jet_iter!= PFJets_corr->end(); ++jet_iter){
181 <      
182 <      reco::PFJet myJet = reco::PFJet(*jet_iter);
183 <      
184 <      if(myJet.et() > 5){
185 <        
186 <        new ((*jetP4_ak5PF_corr)[PFcount_corr]) TLorentzVector(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
187 <        new ((*jetVtx3_ak5PF_corr)[PFcount_corr]) TVector3(myJet.vx(), myJet.vy(), myJet.vz());
188 <        
189 <        ++PFcount_corr;
190 <      }
191 <    }
192 <  }
193 <  
194 <  if(doCaloJets_){
195 <    
196 <    Calocount = 0;
197 <    
198 <    Handle<reco::CaloJetCollection> CaloJets;
199 <    iEvent.getByLabel(CaloJetHandle_, CaloJets);
200 <    
201 <    for(CaloJetCollection::const_iterator jet_iter = CaloJets->begin(); jet_iter!= CaloJets->end(); ++jet_iter){
202 <      
203 <      reco::CaloJet myJet = reco::CaloJet(*jet_iter);
204 <      
174 >
175 >      float scale2 = correctorL2->correction(myJet);
176 >      float scale3 = correctorL3->correction(myJet);
177 >
178        if(myJet.et() > 5){
179          
180 <        new ((*jetP4_ak5Calo)[Calocount]) TLorentzVector(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
208 <        new ((*jetVtx3_ak5Calo)[Calocount]) TVector3(myJet.vx(), myJet.vy(), myJet.vz());
209 <        
210 <        ++Calocount;
180 >        TCJet jetCon;
181          
182 <      }
183 <    }
184 <    
185 <    Handle<reco::CaloJetCollection> CaloJets_corr;
186 <    iEvent.getByLabel(CaloJetHandle_corr, CaloJets_corr);
187 <    
188 <    Calocount_corr = 0;
189 <    
190 <    for(CaloJetCollection::const_iterator jet_iter = CaloJets_corr->begin(); jet_iter!= CaloJets_corr->end(); ++jet_iter){
191 <      
192 <      reco::CaloJet myJet = reco::CaloJet(*jet_iter);
193 <      
194 <      if(myJet.et() > 5){
195 <        
196 <        new ((*jetP4_ak5Calo_corr)[Calocount_corr]) TLorentzVector(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
227 <        new ((*jetVtx3_ak5Calo_corr)[Calocount_corr]) TVector3(myJet.vx(), myJet.vy(), myJet.vz());
182 >        jetCon.SetP4(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
183 >        jetCon.SetVtx(myJet.vx(), myJet.vy(), myJet.vz());
184 >
185 >        jetCon.SetChHadFrac(myJet.chargedHadronEnergyFraction());
186 >        jetCon.SetNeuHadFrac(myJet.neutralHadronEnergyFraction());
187 >        jetCon.SetChEmFrac(myJet.chargedEmEnergyFraction());
188 >        jetCon.SetNeuEmFrac(myJet.neutralEmEnergyFraction());
189 >
190 >        jetCon.SetNumConstit(myJet.chargedMultiplicity() + myJet.neutralMultiplicity());
191 >        jetCon.SetNumChPart(myJet.chargedMultiplicity());
192 >
193 >        jetCon.SetNumChPart(myJet.chargedMultiplicity());
194 >
195 >        jetCon.SetJetCorr(2, scale2);
196 >        jetCon.SetJetCorr(3, scale3);
197          
198 <        ++Calocount_corr;
199 <      }
200 <    }
198 >        new ((*jet_ak5PF)[PFcount]) TCJet(jetCon);
199 >                
200 >        ++PFcount;
201 >      }      
202 >    }  
203    }
204    
205    if(doGenJets_){
# Line 241 | Line 212 | void MPIntuple::analyze(const edm::Event
212        
213        reco::GenJet myJet = reco::GenJet(*jet_iter);
214        
215 <      if(myJet.et() > 5){
215 >      if(myJet.pt() > 5){
216          
217          new ((*jetP4_ak5Gen)[Gencount]) TLorentzVector(myJet.px(), myJet.py(), myJet.pz(), myJet.energy());
247        new ((*jetVtx3_ak5Gen)[Gencount]) TVector3(myJet.vx(), myJet.vy(), myJet.vz());
218          
219          ++Gencount;
220          
221        }
222      }
223    }
224 <  
225 <  
224 >
225 >
226    Handle<reco::VertexCollection> primaryVtx;
227    iEvent.getByLabel(PrimaryVtxHandle_, primaryVtx);
228  
# Line 260 | Line 230 | void MPIntuple::analyze(const edm::Event
230      
231      reco::Vertex myVtx = reco::Vertex(*vtx_iter);
232      
233 <    new ((*PrimaryVtx3)[Vtxcount]) TVector3(myVtx.x(), myVtx.y(), myVtx.z());
234 <    
233 >    TCPrimaryVtx vtxCon;
234 >      
235 >    vtxCon.SetPosition(myVtx.x(), myVtx.y(), myVtx.z());
236 >    vtxCon.SetNDof(myVtx.ndof());
237 >    vtxCon.SetChi2(myVtx.chi2());
238 >      
239 >    new ((*vtx_offline)[Vtxcount]) TCPrimaryVtx(vtxCon);
240 >      
241      ++Vtxcount;
242      
243    }
244 +
245 +
246 +
247 +  //----------  Filling HLT trigger bits!  ------------
248 +
249 +  edm::Handle<TriggerResults> hltR;
250 +  triggerResultsTag_ = InputTag(hlTriggerResults_,"",hltName_);
251 +  iEvent.getByLabel(triggerResultsTag_,hltR);
252 +
253 +  const TriggerNames & triggerNames = iEvent.triggerNames(*hltR);
254 +  hlNames=triggerNames.triggerNames();
255 +
256 +  string MPI_TriggerNames[32] = {"HLT_PixelTracks_Multiplicity70", "HLT_MinBiasBSC_NoBPTX", "HLT_PixelTracks_Multiplicity40","HLT_L1Tech_HCAL_HF", "HLT_IsoTrackHB_8E29", "HLT_IsoTrackHE_8E29", "HLT_L1Tech_RPC_TTU_RBst1_collisions", "HLT_L1_BscMinBiasOR_BptxPlusORMinus", "HLT_L1Tech_BSC_halo_forPhysicsBackground", "HLT_L1Tech_BSC_HighMultiplicity", "HLT_MinBiasPixel_DoubleIsoTrack5", "HLT_MinBiasPixel_DoubleTrack", "HLT_MinBiasPixel_SingleTrack", "HLT_ZeroBiasPixel_SingleTrack", "HLT_MinBiasBSC", "HLT_StoppedHSCP_8E29", "HLT_Jet15U_HcalNoiseFiltered", "HLT_QuadJet15U", "HLT_DiJetAve30U_8E29", "HLT_DiJetAve15U_8E29", "HLT_FwdJet20U", "HLT_Jet50U", "HLT_Jet30U", "HLT_Jet15U", "HLT_BTagMu_Jet10U", "HLT_DoubleJet15U_ForwardBackward", "HLT_BTagIP_Jet50U", "HLT_DoubleLooseIsoTau15", "HLT_SingleLooseIsoTau20", "HLT_HT100U", "HLT_MET100", "HLT_MET45"};
257 +
258 +  bool triggerPassed = false;
259 +  unsigned int triggerStatus = 0x0;
260    
261 +  for (uint iT=0; iT<hlNames.size(); ++iT) {
262 +
263 +    triggerPassed = triggerDecision(hltR, iT);
264 +    
265 +    if(triggerPassed){
266 +      
267 +      for (int j = 0; j != 32; ++j){
268 +        
269 +        if (hlNames[iT] == MPI_TriggerNames[j])
270 +          {
271 +            cout<<"trigger name: "<<hlNames[iT]<<"  status: "<<triggerPassed<<endl;
272 +            triggerStatus |= 0x01 << j;
273 +          }
274 +      }
275 +    }
276 +  }
277    
278 <  if(Calocount > 0 || Calocount_corr >0 || PFcount_corr > 0 || PFcount > 0 || Gencount > 0)  sTree -> Fill();
278 >  cout<<triggerStatus<<endl;
279  
280 <  jetP4_ak5PF->Clear();
273 <  jetP4_ak5PF_corr->Clear();
274 <  jetVtx3_ak5PF->Clear();
275 <  jetP4_ak5Calo->Clear();
276 <  jetP4_ak5Calo_corr->Clear();
277 <  jetVtx3_ak5Calo->Clear();
278 <  jetP4_ak5Gen->Clear();
279 <  jetVtx3_ak5Gen->Clear();
280 >  if(PFcount > 0 || Gencount > 0 || Vtxcount > 0);  sTree -> Fill();
281    
282 +  jet_ak5PF->Clear();
283 +  vtx_offline->Clear();
284 +  jetP4_ak5Gen->Clear();
285 +
286   }
287  
288  
# Line 288 | Line 293 | void  MPIntuple::beginJob()
293    ntupleFile           = new TFile(rootfilename.c_str(), "RECREATE");
294    sTree                = new TTree("dpsTree", "Tree for Jets");
295    
296 <  jetP4_ak5PF          = new TClonesArray("TLorentzVector");
292 <  jetVtx3_ak5PF        = new TClonesArray("TVector3");
293 <  jetP4_ak5Calo        = new TClonesArray("TLorentzVector");
294 <  jetVtx3_ak5Calo      = new TClonesArray("TVector3");
295 <  jetP4_ak5PF_corr     = new TClonesArray("TLorentzVector");
296 <  jetVtx3_ak5PF_corr   = new TClonesArray("TVector3");
297 <  jetP4_ak5Calo_corr   = new TClonesArray("TLorentzVector");
298 <  jetVtx3_ak5Calo_corr = new TClonesArray("TVector3");
296 >  jet_ak5PF            = new TClonesArray("TCJet");
297    jetP4_ak5Gen         = new TClonesArray("TLorentzVector");
298 <  jetVtx3_ak5Gen       = new TClonesArray("TVector3");
301 <  PrimaryVtx3          = new TClonesArray("TVector3");
298 >  vtx_offline          = new TClonesArray("TCPrimaryVtx");
299  
300 <  sTree->Branch("jetP4_ak5PF",&jetP4_ak5PF, 6400, 0);
304 <  sTree->Branch("jetP4_ak5PF_corr",&jetP4_ak5PF_corr, 6400, 0);
305 <  sTree->Branch("jetVtx3_ak5PF",&jetVtx3_ak5PF, 6400, 0);
306 <  sTree->Branch("jetVtx3_ak5PF_corr",&jetVtx3_ak5PF, 6400, 0);
307 <  sTree->Branch("jetP4_ak5Calo",&jetP4_ak5Calo, 6400, 0);
308 <  sTree->Branch("jetP4_ak5Calo_corr",&jetP4_ak5Calo_corr, 6400, 0);
309 <  sTree->Branch("jetVtx3_ak5Calo",&jetVtx3_ak5Calo, 6400, 0);
310 <  sTree->Branch("jetVtx3_ak5Calo_corr",&jetVtx3_ak5Calo, 6400, 0);
300 >  sTree->Branch("jet_ak5PF",&jet_ak5PF, 6400, 0);
301    sTree->Branch("jetP4_ak5Gen",&jetP4_ak5Gen, 6400, 0);
302 <  sTree->Branch("jetVtx3_ak5Gen",&jetVtx3_ak5Gen, 6400, 0);
313 <  sTree->Branch("PrimaryVtx",&PrimaryVtx3, 6400, 0);
302 >  sTree->Branch("vtx_offline",&vtx_offline, 6400, 0);
303  
304    sTree->Branch("eventNumber",&eventNumber, "eventNumber/I");
305    sTree->Branch("runNumber",&runNumber, "runNumber/I");
# Line 327 | Line 316 | void MPIntuple::endJob()
316   }
317  
318  
319 + bool MPIntuple::triggerDecision(edm::Handle<edm::TriggerResults> &hltR, int iTrigger)
320 + {
321 +    bool triggerPassed = false;
322 +    if(hltR->wasrun(iTrigger) &&
323 +       hltR->accept(iTrigger) &&
324 +       !hltR->error(iTrigger) ){
325 +      triggerPassed = true;
326 +    }
327 +    return triggerPassed;
328 +  }
329 +
330 +
331   //define this as a plug-in
332   DEFINE_FWK_MODULE(MPIntuple);

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines