ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MuJetAnalysis/AnalysisTools/IsolationStudyNtuple/src/IsolationStudyNtuple.cc
Revision: 1.2
Committed: Wed Nov 10 01:16:03 2010 UTC (14 years, 5 months ago) by pivarski
Content type: text/plain
Branch: MAIN
Changes since 1.1: +213 -15 lines
Log Message:
add a minimum threshold to track isolation

File Contents

# User Rev Content
1 pivarski 1.1 // -*- C++ -*-
2     //
3     // Package: IsolationStudyNtuple
4     // Class: IsolationStudyNtuple
5     //
6     /**\class IsolationStudyNtuple IsolationStudyNtuple.cc AnalysisTools/IsolationStudyNtuple/src/IsolationStudyNtuple.cc
7    
8     Description: [one line class summary]
9    
10     Implementation:
11     [Notes on implementation]
12     */
13     //
14     // Original Author: James Pivarski
15     // Created: Tue Nov 9 13:34:12 CST 2010
16 pivarski 1.2 // $Id: IsolationStudyNtuple.cc,v 1.1 2010/11/09 19:37:57 pivarski Exp $
17 pivarski 1.1 //
18     //
19    
20    
21     // system include files
22     #include <memory>
23    
24     // user include files
25     #include "FWCore/Framework/interface/Frameworkfwd.h"
26     #include "FWCore/Framework/interface/EDAnalyzer.h"
27     #include "FWCore/Framework/interface/Event.h"
28     #include "FWCore/Framework/interface/MakerMacros.h"
29 pivarski 1.2 #include "FWCore/ParameterSet/interface/ParameterSet.h"
30    
31     #include "AnalysisDataFormats/MuJetAnalysis/interface/MultiMuon.h"
32     #include "CommonTools/UtilAlgos/interface/TFileService.h"
33     #include "DataFormats/MuonReco/interface/Muon.h"
34     #include "DataFormats/PatCandidates/interface/Muon.h"
35     #include "DataFormats/TrackReco/interface/TrackFwd.h"
36     #include "DataFormats/PatCandidates/interface/TriggerEvent.h"
37     #include "DataFormats/PatCandidates/interface/TriggerPath.h"
38     #include "FWCore/ServiceRegistry/interface/Service.h"
39     #include "TTree.h"
40 pivarski 1.1
41     //
42     // class declaration
43     //
44    
45     class IsolationStudyNtuple : public edm::EDAnalyzer {
46     public:
47     explicit IsolationStudyNtuple(const edm::ParameterSet&);
48     ~IsolationStudyNtuple();
49    
50    
51     private:
52     virtual void beginJob() ;
53     virtual void analyze(const edm::Event&, const edm::EventSetup&);
54     virtual void endJob() ;
55    
56     // ----------member data ---------------------------
57 pivarski 1.2 TTree *m_iso01th00, *m_iso03th00, *m_iso05th00, *m_iso07th00, *m_iso10th00;
58     TTree *m_iso01th01, *m_iso03th01, *m_iso05th01, *m_iso07th01, *m_iso10th01;
59     TTree *m_iso01th02, *m_iso03th02, *m_iso05th02, *m_iso07th02, *m_iso10th02;
60     TTree *m_iso01th05, *m_iso03th05, *m_iso05th05, *m_iso07th05, *m_iso10th05;
61     TTree *m_iso01th10, *m_iso03th10, *m_iso05th10, *m_iso07th10, *m_iso10th10;
62     TTree *m_iso01th20, *m_iso03th20, *m_iso05th20, *m_iso07th20, *m_iso10th20;
63     Int_t m_daughters;
64     Float_t m_mass;
65     Float_t m_dRmax;
66     Float_t m_dphimax;
67     Float_t m_detamax;
68     Float_t m_vertexProb;
69     Float_t m_pt;
70     Float_t m_phi;
71     Float_t m_eta;
72     Float_t m_pt1;
73     Float_t m_pt2;
74     Float_t m_trackiso;
75     Int_t m_numAboveThresh;
76 pivarski 1.1 };
77    
78     //
79     // constants, enums and typedefs
80     //
81    
82     //
83     // static data member definitions
84     //
85    
86     //
87     // constructors and destructor
88     //
89     IsolationStudyNtuple::IsolationStudyNtuple(const edm::ParameterSet& iConfig)
90    
91     {
92     //now do what ever initialization is needed
93 pivarski 1.2 edm::Service<TFileService> tFile;
94     m_iso01th00 = tFile->make<TTree>("iso01th00", "iso01th00"); m_iso03th00 = tFile->make<TTree>("iso03th00", "iso03th00"); m_iso05th00 = tFile->make<TTree>("iso05th00", "iso05th00"); m_iso07th00 = tFile->make<TTree>("iso07th00", "iso07th00"); m_iso10th00 = tFile->make<TTree>("iso10th00", "iso10th00");
95     m_iso01th01 = tFile->make<TTree>("iso01th01", "iso01th01"); m_iso03th01 = tFile->make<TTree>("iso03th01", "iso03th01"); m_iso05th01 = tFile->make<TTree>("iso05th01", "iso05th01"); m_iso07th01 = tFile->make<TTree>("iso07th01", "iso07th01"); m_iso10th01 = tFile->make<TTree>("iso10th01", "iso10th01");
96     m_iso01th02 = tFile->make<TTree>("iso01th02", "iso01th02"); m_iso03th02 = tFile->make<TTree>("iso03th02", "iso03th02"); m_iso05th02 = tFile->make<TTree>("iso05th02", "iso05th02"); m_iso07th02 = tFile->make<TTree>("iso07th02", "iso07th02"); m_iso10th02 = tFile->make<TTree>("iso10th02", "iso10th02");
97     m_iso01th05 = tFile->make<TTree>("iso01th05", "iso01th05"); m_iso03th05 = tFile->make<TTree>("iso03th05", "iso03th05"); m_iso05th05 = tFile->make<TTree>("iso05th05", "iso05th05"); m_iso07th05 = tFile->make<TTree>("iso07th05", "iso07th05"); m_iso10th05 = tFile->make<TTree>("iso10th05", "iso10th05");
98     m_iso01th10 = tFile->make<TTree>("iso01th10", "iso01th10"); m_iso03th10 = tFile->make<TTree>("iso03th10", "iso03th10"); m_iso05th10 = tFile->make<TTree>("iso05th10", "iso05th10"); m_iso07th10 = tFile->make<TTree>("iso07th10", "iso07th10"); m_iso10th10 = tFile->make<TTree>("iso10th10", "iso10th10");
99     m_iso01th20 = tFile->make<TTree>("iso01th20", "iso01th20"); m_iso03th20 = tFile->make<TTree>("iso03th20", "iso03th20"); m_iso05th20 = tFile->make<TTree>("iso05th20", "iso05th20"); m_iso07th20 = tFile->make<TTree>("iso07th20", "iso07th20"); m_iso10th20 = tFile->make<TTree>("iso10th20", "iso10th20");
100    
101     for (int cone = 1; cone <= 5; cone++) {
102     for (int thresh = 0; thresh <= 5; thresh++) {
103     TTree *ttree = NULL;
104     if (cone == 1) {
105     if (thresh == 0) ttree = m_iso01th00;
106     else if (thresh == 1) ttree = m_iso01th01;
107     else if (thresh == 2) ttree = m_iso01th02;
108     else if (thresh == 3) ttree = m_iso01th05;
109     else if (thresh == 4) ttree = m_iso01th10;
110     else if (thresh == 5) ttree = m_iso01th20;
111     else assert(false);
112     }
113     else if (cone == 2) {
114     if (thresh == 0) ttree = m_iso03th00;
115     else if (thresh == 1) ttree = m_iso03th01;
116     else if (thresh == 2) ttree = m_iso03th02;
117     else if (thresh == 3) ttree = m_iso03th05;
118     else if (thresh == 4) ttree = m_iso03th10;
119     else if (thresh == 5) ttree = m_iso03th20;
120     else assert(false);
121     }
122     else if (cone == 3) {
123     if (thresh == 0) ttree = m_iso05th00;
124     else if (thresh == 1) ttree = m_iso05th01;
125     else if (thresh == 2) ttree = m_iso05th02;
126     else if (thresh == 3) ttree = m_iso05th05;
127     else if (thresh == 4) ttree = m_iso05th10;
128     else if (thresh == 5) ttree = m_iso05th20;
129     else assert(false);
130     }
131     else if (cone == 4) {
132     if (thresh == 0) ttree = m_iso07th00;
133     else if (thresh == 1) ttree = m_iso07th01;
134     else if (thresh == 2) ttree = m_iso07th02;
135     else if (thresh == 3) ttree = m_iso07th05;
136     else if (thresh == 4) ttree = m_iso07th10;
137     else if (thresh == 5) ttree = m_iso07th20;
138     else assert(false);
139     }
140     else if (cone == 5) {
141     if (thresh == 0) ttree = m_iso10th00;
142     else if (thresh == 1) ttree = m_iso10th01;
143     else if (thresh == 2) ttree = m_iso10th02;
144     else if (thresh == 3) ttree = m_iso10th05;
145     else if (thresh == 4) ttree = m_iso10th10;
146     else if (thresh == 5) ttree = m_iso10th20;
147     else assert(false);
148     }
149     else assert(false);
150    
151     ttree->Branch("daughters", &m_daughters, "daughters/I");
152     ttree->Branch("mass", &m_mass, "mass/F");
153     ttree->Branch("dRmax", &m_dRmax, "dRmax/F");
154     ttree->Branch("dphimax", &m_dphimax, "dphimax/F");
155     ttree->Branch("detamax", &m_detamax, "detamax/F");
156     ttree->Branch("vertexProb", &m_vertexProb, "vertexProb/F");
157     ttree->Branch("pt", &m_pt, "pt/F");
158     ttree->Branch("phi", &m_phi, "phi/F");
159     ttree->Branch("eta", &m_eta, "eta/F");
160     ttree->Branch("pt1", &m_pt1, "pt1/F");
161     ttree->Branch("pt2", &m_pt2, "pt2/F");
162     ttree->Branch("trackiso", &m_trackiso, "trackiso/F");
163     ttree->Branch("numAboveThresh", &m_numAboveThresh, "numAboveThresh/I");
164     }
165     }
166 pivarski 1.1 }
167    
168    
169     IsolationStudyNtuple::~IsolationStudyNtuple()
170     {
171    
172     // do anything here that needs to be done at desctruction time
173     // (e.g. close files, deallocate resources etc.)
174    
175     }
176    
177    
178     //
179     // member functions
180     //
181    
182     // ------------ method called to for each event ------------
183     void
184     IsolationStudyNtuple::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
185     {
186 pivarski 1.2 edm::Handle<pat::MultiMuonCollection> muJets;
187     iEvent.getByLabel("MuJetProducer", muJets);
188 pivarski 1.1
189 pivarski 1.2 edm::Handle<reco::TrackCollection> tracks;
190     iEvent.getByLabel("generalTracks", tracks);
191 pivarski 1.1
192 pivarski 1.2 for (int cone = 1; cone <= 5; cone++) {
193     for (int thresh = 0; thresh <= 5; thresh++) {
194     double centralCone;
195     if (cone == 1) centralCone = 0.1;
196     else if (cone == 2) centralCone = 0.3;
197     else if (cone == 3) centralCone = 0.5;
198     else if (cone == 4) centralCone = 0.7;
199     else if (cone == 5) centralCone = 1.0;
200     else assert(false);
201    
202     double centralThreshold;
203     if (thresh == 0) centralThreshold = 0.0;
204     else if (thresh == 1) centralThreshold = 0.1;
205     else if (thresh == 2) centralThreshold = 0.2;
206     else if (thresh == 3) centralThreshold = 0.5;
207     else if (thresh == 4) centralThreshold = 1.0;
208     else if (thresh == 5) centralThreshold = 2.0;
209     else assert(false);
210    
211     TTree *ttree = NULL;
212     if (cone == 1) {
213     if (thresh == 0) ttree = m_iso01th00;
214     else if (thresh == 1) ttree = m_iso01th01;
215     else if (thresh == 2) ttree = m_iso01th02;
216     else if (thresh == 3) ttree = m_iso01th05;
217     else if (thresh == 4) ttree = m_iso01th10;
218     else if (thresh == 5) ttree = m_iso01th20;
219     else assert(false);
220     }
221     else if (cone == 2) {
222     if (thresh == 0) ttree = m_iso03th00;
223     else if (thresh == 1) ttree = m_iso03th01;
224     else if (thresh == 2) ttree = m_iso03th02;
225     else if (thresh == 3) ttree = m_iso03th05;
226     else if (thresh == 4) ttree = m_iso03th10;
227     else if (thresh == 5) ttree = m_iso03th20;
228     else assert(false);
229     }
230     else if (cone == 3) {
231     if (thresh == 0) ttree = m_iso05th00;
232     else if (thresh == 1) ttree = m_iso05th01;
233     else if (thresh == 2) ttree = m_iso05th02;
234     else if (thresh == 3) ttree = m_iso05th05;
235     else if (thresh == 4) ttree = m_iso05th10;
236     else if (thresh == 5) ttree = m_iso05th20;
237     else assert(false);
238     }
239     else if (cone == 4) {
240     if (thresh == 0) ttree = m_iso07th00;
241     else if (thresh == 1) ttree = m_iso07th01;
242     else if (thresh == 2) ttree = m_iso07th02;
243     else if (thresh == 3) ttree = m_iso07th05;
244     else if (thresh == 4) ttree = m_iso07th10;
245     else if (thresh == 5) ttree = m_iso07th20;
246     else assert(false);
247     }
248     else if (cone == 5) {
249     if (thresh == 0) ttree = m_iso10th00;
250     else if (thresh == 1) ttree = m_iso10th01;
251     else if (thresh == 2) ttree = m_iso10th02;
252     else if (thresh == 3) ttree = m_iso10th05;
253     else if (thresh == 4) ttree = m_iso10th10;
254     else if (thresh == 5) ttree = m_iso10th20;
255     else assert(false);
256     }
257     else assert(false);
258    
259     for (pat::MultiMuonCollection::const_iterator muJetOriginal = muJets->begin(); muJetOriginal != muJets->end(); ++muJetOriginal) {
260     pat::MultiMuon muJet(*muJetOriginal);
261     muJet.calculateTrackIsolation(&*tracks, centralCone, 0., centralThreshold, 0.);
262     muJet.calculateNumberAboveThresholdIsolation(&*tracks, centralCone, 0., centralThreshold, 0.);
263    
264     m_daughters = muJet.numberOfDaughters();
265     m_mass = muJet.mass();
266     m_dRmax = muJet.dRmax();
267     m_dphimax = 0.;
268     m_detamax = 0.;
269     m_pt1 = 0.;
270     m_pt2 = 0.;
271     for (int i = 0; i < m_daughters; i++) {
272     for (int j = i+1; j < m_daughters; j++) {
273     if (fabs(muJet.dphi(i, j)) > fabs(m_dphimax)) m_dphimax = muJet.dphi(i, j);
274     if (fabs(muJet.deta(i, j)) > fabs(m_detamax)) m_detamax = muJet.deta(i, j);
275     }
276    
277     if (muJet.daughter(i)->pt() > m_pt1) {
278     m_pt2 = m_pt1;
279     m_pt1 = muJet.daughter(i)->pt();
280     }
281     else if (muJet.daughter(i)->pt() > m_pt2) {
282     m_pt2 = muJet.daughter(i)->pt();
283     }
284     }
285     if (muJet.vertexValid()) m_vertexProb = muJet.vertexProb();
286     else m_vertexProb = -1.;
287     m_pt = muJet.pt();
288     m_phi = muJet.phi();
289     m_eta = muJet.eta();
290     m_trackiso = muJet.centralTrackIsolation();
291     m_numAboveThresh = muJet.centralNumberAboveThreshold();
292    
293     ttree->Fill();
294     }
295     }
296     }
297 pivarski 1.1 }
298    
299    
300     // ------------ method called once each job just before starting event loop ------------
301     void
302     IsolationStudyNtuple::beginJob()
303     {
304     }
305    
306     // ------------ method called once each job just after ending the event loop ------------
307     void
308     IsolationStudyNtuple::endJob() {
309     }
310    
311     //define this as a plug-in
312     DEFINE_FWK_MODULE(IsolationStudyNtuple);