ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/CmsHi/JetAnalysis/src/HiPFJetAnalyzer.cc
Revision: 1.1
Committed: Tue Sep 20 19:45:06 2011 UTC (13 years, 7 months ago) by yilmaz
Content type: text/plain
Branch: MAIN
CVS Tags: HiForest_V02_85, HiForest_V02_84, HiForest_V02_83, HiForest_V02_82, HiForest_V02_81, HiForest_V02_80, HiForest_V02_79, HiForest_V02_78, HiForest_V02_77, HiForest_V02_76, HiForest_V02_75, HiForest_V02_74, HiForest_V02_73, HiForest_V02_72, HiForest_V02_71, HiForest_V02_70, HiForest_V02_69, HiForest_V02_68, HiForest_V02_67, HiForest_V02_66, HiForest_V02_65, HiForest_V02_64, HiForest_V02_63, HiForest_V02_62, HiForest_V02_61, HiForest_V02_60, HiForest_V02_59, HiForest_V02_58, HiForest_V02_57, HiForest_V02_56, HiForest_V02_55, HiForest_V02_54, HiForest_V02_53, HiForest_V02_52, HiForest_V02_51, HiForest_V02_50, HiForest_V02_49, HiForest_V02_48, HiForest_V02_47, HiForest_V02_46, HiForest_V02_45, HiForest_V02_44, HiForest_V02_43, HiForest_V02_42, HiForest_V02_41, HiForest_V02_40, HiForest_V02_39, HiForest_V02_38, HiForest_V02_37, HiForest_V02_36, HiForest_V02_35, HiForest_V02_34, HiForest_V02_33, HiForest_V02_32, HiForest_V02_31, HiForest_V02_30, HiForest_V02_27, HiForest_V02_26, QM_2012, HiForest_V02_25, HiForest_V02_24, HiForest_V02_23, HiForest_V02_22, HiForest_V02_21, HiForest_V02_20, HiForest_V02_19, HiForest_V02_18, HiForest_V02_17, HiForest_V02_16, HiForest_V02_15, HiForest_V02_14, HiForest_V02_13, HiForest_V02_12, HiForest_V02_11, HiForest_V02_10, HiForest_V02_09, HiForest_V02_08, HiForest_V02_07, HiForest_V02_06, HiForest_V02_05, HiForest_V02_04, HiForest_V02_03, HiForest_V02_02, HiForest_V02_01, HiForest_V02_00, hi44X_02, hi413_03, hi413_02, hi441_1, hi441_0, hi413_11, hi413_10, hi413_09, HEAD
Branch point for: branch_44x
Log Message:
combining all analyzers

File Contents

# User Rev Content
1 yilmaz 1.1 /*
2     Based on the jet response analyzer
3     Modified by Matt Nguyen, November 2010
4    
5     */
6    
7     #include "CmsHi/JetAnalysis/interface/HiPFJetAnalyzer.h"
8    
9     #include "FWCore/Framework/interface/EventSetup.h"
10     #include "FWCore/Framework/interface/ESHandle.h"
11     //#include <Math/DistFunc.h>
12     #include "TMath.h"
13    
14     #include "FWCore/Framework/interface/ESHandle.h"
15    
16     #include "FWCore/MessageLogger/interface/MessageLogger.h"
17     #include "FWCore/Utilities/interface/Exception.h"
18     #include "FWCore/Framework/interface/EventSetup.h"
19    
20    
21     #include "DataFormats/HeavyIonEvent/interface/CentralityBins.h"
22    
23     #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
24     #include "DataFormats/HeavyIonEvent/interface/Centrality.h"
25     #include "DataFormats/PatCandidates/interface/Jet.h"
26     #include "DataFormats/JetReco/interface/CaloJetCollection.h"
27     #include "SimDataFormats/HiGenData/interface/GenHIEvent.h"
28     #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
29     #include "DataFormats/JetReco/interface/GenJetCollection.h"
30     #include "DataFormats/Math/interface/deltaR.h"
31    
32     #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
33     #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
34     #include "SimTracker/Records/interface/TrackAssociatorRecord.h"
35     #include "DataFormats/RecoCandidate/interface/TrackAssociation.h"
36     #include "SimTracker/TrackAssociation/interface/TrackAssociatorByHits.h"
37    
38     #include "DataFormats/JetReco/interface/CaloJetCollection.h"
39     #include "DataFormats/JetReco/interface/PFJetCollection.h"
40     #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
41     #include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
42     #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
43     #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
44     #include "DataFormats/TrackReco/interface/Track.h"
45     #include "DataFormats/TrackReco/interface/TrackFwd.h"
46    
47     #include "DataFormats/Common/interface/TriggerResults.h"
48     #include "FWCore/Common/interface/TriggerNames.h"
49     #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetup.h"
50     #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
51     #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMapRecord.h"
52     #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerObjectMap.h"
53     #include "L1Trigger/GlobalTrigger/interface/L1GlobalTrigger.h"
54    
55    
56     #include "DataFormats/BeamSpot/interface/BeamSpot.h"
57    
58    
59     using namespace std;
60     using namespace edm;
61     using namespace reco;
62    
63    
64     HiPFJetAnalyzer::HiPFJetAnalyzer(const edm::ParameterSet& iConfig) {
65    
66    
67     jetTag1_ = iConfig.getParameter<InputTag>("jetTag1");
68     jetTag2_ = iConfig.getParameter<InputTag>("jetTag2");
69     jetTag3_ = iConfig.getParameter<InputTag>("jetTag3");
70     jetTag4_ = iConfig.getParameter<InputTag>("jetTag4");
71    
72     recoJetTag1_ = iConfig.getParameter<InputTag>("recoJetTag1");
73     recoJetTag2_ = iConfig.getParameter<InputTag>("recoJetTag2");
74     recoJetTag3_ = iConfig.getParameter<InputTag>("recoJetTag3");
75     recoJetTag4_ = iConfig.getParameter<InputTag>("recoJetTag4");
76    
77     genJetTag1_ = iConfig.getParameter<InputTag>("genJetTag1");
78     genJetTag2_ = iConfig.getParameter<InputTag>("genJetTag2");
79     genJetTag3_ = iConfig.getParameter<InputTag>("genJetTag3");
80     genJetTag4_ = iConfig.getParameter<InputTag>("genJetTag4");
81    
82     pfCandidatesTag_ = iConfig.getParameter<InputTag>("pfCandidatesTag");
83     trackTag_ = iConfig.getParameter<edm::InputTag>("trackTag");
84     vertexTag_ = iConfig.getParameter<edm::InputTag>("vertexTag");
85    
86     verbose_ = iConfig.getUntrackedParameter<bool>("verbose",false);
87    
88     useCentrality_ = iConfig.getUntrackedParameter<bool>("useCentrality",false);
89     isMC_ = iConfig.getUntrackedParameter<bool>("isMC",false);
90     genParticleThresh_ = iConfig.getParameter<double>("genParticleThresh");
91    
92     genParticleTag_ = iConfig.getParameter<InputTag>("genParticleTag");
93     eventInfoTag_ = iConfig.getParameter<InputTag>("eventInfoTag");
94    
95     hasSimInfo_ = iConfig.getUntrackedParameter<bool>("hasSimInfo");
96     simTracksTag_ = iConfig.getParameter<InputTag>("SimTracks");
97    
98     if(!isMC_){
99     L1gtReadout_ = iConfig.getParameter<edm::InputTag>("L1gtReadout");
100     hltResName_ = iConfig.getUntrackedParameter<string>("hltTrgResults","TriggerResults::HLT");
101    
102    
103     if (iConfig.exists("hltTrgNames"))
104     hltTrgNames_ = iConfig.getUntrackedParameter<vector<string> >("hltTrgNames");
105    
106     if (iConfig.exists("hltProcNames"))
107     hltProcNames_ = iConfig.getUntrackedParameter<vector<string> >("hltProcNames");
108     else {
109     hltProcNames_.push_back("FU");
110     hltProcNames_.push_back("HLT");
111     }
112     }
113    
114    
115     cout<<" tracks : "<<trackTag_<<endl;
116     cout<<" jet collection 1: "<<jetTag1_<<endl;
117     cout<<" jet collection 2: "<<jetTag2_<<endl;
118     cout<<" jet collection 3: "<<jetTag3_<<endl;
119     cout<<" jet collection 4: "<<jetTag4_<<endl;
120    
121     jets_.nj1 = 0;
122     jets_.nj2 = 0;
123     jets_.nj3 = 0;
124     jets_.nj4 = 0;
125     jets_.nunmatch_j1 = 0;
126     jets_.nunmatch_j2 = 0;
127     jets_.nunmatch_j3 = 0;
128     jets_.nunmatch_j4 = 0;
129     jets_.nPFcand = 0;
130     jets_.ntrack = 0;
131     jets_.ngenp = 0;
132    
133    
134    
135     }
136    
137    
138    
139     HiPFJetAnalyzer::~HiPFJetAnalyzer() { }
140    
141    
142    
143     void
144     HiPFJetAnalyzer::beginRun(const edm::Run& run,
145     const edm::EventSetup & es) {}
146    
147     void
148     HiPFJetAnalyzer::beginJob() {
149    
150     centrality_ = 0;
151    
152    
153     t = fs2->make<TTree>("t","Jet Analysis Tree");
154    
155    
156     // TTree* t= new TTree("t","Jet Response Analyzer");
157     t->Branch("run",&jets_.run,"run/I");
158     t->Branch("evt",&jets_.evt,"evt/I");
159     t->Branch("lumi",&jets_.lumi,"lumi/I");
160     t->Branch("b",&jets_.b,"b/F");
161    
162     t->Branch("vx",&jets_.vx,"vx/F");
163     t->Branch("vy",&jets_.vy,"vy/F");
164     t->Branch("vz",&jets_.vz,"vz/F");
165    
166     t->Branch("hasVtx",&jets_.hasVtx,"hasVtx/O");
167    
168     t->Branch("vxErr",&jets_.vxErr,"vxErr/F");
169     t->Branch("vyErr",&jets_.vyErr,"vyErr/F");
170     t->Branch("vzErr",&jets_.vzErr,"vzErr/F");
171    
172     t->Branch("bsx",&jets_.bsx,"bsx/F");
173     t->Branch("bsy",&jets_.bsy,"bsy/F");
174     t->Branch("bsz",&jets_.bsz,"bsz/F");
175     t->Branch("bswx",&jets_.bswx,"bswx/F");
176     t->Branch("bswy",&jets_.bswy,"bswy/F");
177    
178     t->Branch("hf",&jets_.hf,"hf/F");
179     t->Branch("bin",&jets_.bin,"bin/I");
180    
181     // ICPU5
182     t->Branch("nj1",&jets_.nj1,"nj1/I");
183     t->Branch("nj2",&jets_.nj2,"nj2/I");
184     t->Branch("nj3",&jets_.nj3,"nj3/I");
185     t->Branch("nj4",&jets_.nj4,"nj4/I");
186    
187     t->Branch("rawpt_j1",jets_.rawpt_j1,"rawpt_j1[nj1]/F");
188     t->Branch("jtpt_j1",jets_.jtpt_j1,"jtpt_j1[nj1]/F");
189     t->Branch("jteta_j1",jets_.jteta_j1,"jteta_j1[nj1]/F");
190     t->Branch("jty_j1",jets_.jty_j1,"jty_j1[nj1]/F");
191     t->Branch("jtphi_j1",jets_.jtphi_j1,"jtphi_j1[nj1]/F");
192     t->Branch("preL1et_j1",jets_.preL1et_j1,"preL1et_j1[nj1]/F");
193     t->Branch("L2_j1",jets_.L2_j1,"L2_j1[nj1]/F");
194     t->Branch("L3_j1",jets_.L3_j1,"L3_j1[nj1]/F");
195     t->Branch("area_j1",jets_.area_j1,"area_j1[nj1]/F");
196    
197     if(isMC_){
198     // Matched gen jets
199     t->Branch("refpt_j1",jets_.refpt_j1,"refpt_j1[nj1]/F");
200     t->Branch("refeta_j1",jets_.refeta_j1,"refeta_j1[nj1]/F");
201     t->Branch("refy_j1",jets_.refy_j1,"refy_j1[nj1]/F");
202     t->Branch("refphi_j1",jets_.refphi_j1,"refphi_j1[nj1]/F");
203     t->Branch("refdrjt_j1",jets_.refdrjt_j1,"refdrjt_j1[nj1]/F");
204     t->Branch("refpartonpt_j1",jets_.refpartonpt_j1,"refpartonpt_j1[nj1]/F");
205     t->Branch("refpartonflavor_j1",jets_.refpartonflavor_j1,"refpartonflavor_j1[nj1]/F");
206    
207     // Unmatched gen jets
208     t->Branch("nunmatch_j1",&jets_.nunmatch_j1,"nunmatch_j1/I");
209     t->Branch("unmatchpt_j1",jets_.unmatchpt_j1,"unmatchpt_j1[nunmatch_j1]/F");
210     t->Branch("unmatcheta_j1",jets_.unmatcheta_j1,"unmatcheta_j1[nunmatch_j1]/F");
211     t->Branch("unmatchy_j1",jets_.unmatchy_j1,"unmatchy_j1[nunmatch_j1]/F");
212     t->Branch("unmatchphi_j1",jets_.unmatchphi_j1,"unmatchphi_j1[nunmatch_j1]/F");
213     }
214    
215    
216     // J2
217    
218     t->Branch("rawpt_j2",jets_.rawpt_j2,"rawpt_j2[nj2]/F");
219     t->Branch("jtpt_j2",jets_.jtpt_j2,"jtpt_j2[nj2]/F");
220     t->Branch("jteta_j2",jets_.jteta_j2,"jteta_j2[nj2]/F");
221     t->Branch("jty_j2",jets_.jty_j2,"jty_j2[nj2]/F");
222     t->Branch("jtphi_j2",jets_.jtphi_j2,"jtphi_j2[nj2]/F");
223     t->Branch("preL1et_j2",jets_.preL1et_j2,"preL1et_j2[nj2]/F");
224     t->Branch("L2_j2",jets_.L2_j2,"L2_j2[nj2]/F");
225     t->Branch("L3_j2",jets_.L3_j2,"L3_j2[nj2]/F");
226     t->Branch("area_j2",jets_.area_j2,"area_j2[nj2]/F");
227    
228     if(isMC_){
229     t->Branch("refpt_j2",jets_.refpt_j2,"refpt_j2[nj2]/F");
230     t->Branch("refeta_j2",jets_.refeta_j2,"refeta_j2[nj2]/F");
231     t->Branch("refy_j2",jets_.refy_j2,"refy_j2[nj2]/F");
232     t->Branch("refphi_j2",jets_.refphi_j2,"refphi_j2[nj2]/F");
233     t->Branch("refdrjt_j2",jets_.refdrjt_j2,"refdrjt_j2[nj2]/F");
234     t->Branch("refpartonpt_j2",jets_.refpartonpt_j2,"refpartonpt_j2[nj2]/F");
235     t->Branch("refpartonflavor_j2",jets_.refpartonflavor_j2,"refpartonflavor_j2[nj2]/F");
236    
237     // Unmatched gen jets
238     t->Branch("nunmatch_j2",&jets_.nunmatch_j2,"nunmatch_j2/I");
239     t->Branch("unmatchpt_j2",jets_.unmatchpt_j2,"unmatchpt_j2[nunmatch_j2]/F");
240     t->Branch("unmatcheta_j2",jets_.unmatcheta_j2,"unmatcheta_j2[nunmatch_j2]/F");
241     t->Branch("unmatchy_j2",jets_.unmatchy_j2,"unmatchy_j2[nunmatch_j2]/F");
242     t->Branch("unmatchphi_j2",jets_.unmatchphi_j2,"unmatchphi_j2[nunmatch_j2]/F");
243    
244     }
245    
246    
247    
248    
249     // J3
250    
251     t->Branch("rawpt_j3",jets_.rawpt_j3,"rawpt_j3[nj3]/F");
252     t->Branch("jtpt_j3",jets_.jtpt_j3,"jtpt_j3[nj3]/F");
253     t->Branch("jteta_j3",jets_.jteta_j3,"jteta_j3[nj3]/F");
254     t->Branch("jty_j3",jets_.jty_j3,"jty_j3[nj3]/F");
255     t->Branch("jtphi_j3",jets_.jtphi_j3,"jtphi_j3[nj3]/F");
256     t->Branch("preL1et_j3",jets_.preL1et_j3,"preL1et_j3[nj3]/F");
257     t->Branch("L2_j3",jets_.L2_j3,"L2_j3[nj3]/F");
258     t->Branch("L3_j3",jets_.L3_j3,"L3_j3[nj3]/F");
259     t->Branch("area_j3",jets_.area_j3,"area_j3[nj3]/F");
260    
261     if(isMC_){
262     t->Branch("refpt_j3",jets_.refpt_j3,"refpt_j3[nj3]/F");
263     t->Branch("refeta_j3",jets_.refeta_j3,"refeta_j3[nj3]/F");
264     t->Branch("refy_j3",jets_.refy_j3,"refy_j3[nj3]/F");
265     t->Branch("refphi_j3",jets_.refphi_j3,"refphi_j3[nj3]/F");
266     t->Branch("refdrjt_j3",jets_.refdrjt_j3,"refdrjt_j3[nj3]/F");
267     t->Branch("refpartonpt_j3",jets_.refpartonpt_j3,"refpartonpt_j3[nj3]/F");
268     t->Branch("refpartonflavor_j3",jets_.refpartonflavor_j3,"refpartonflavor_j3[nj3]/F");
269    
270     // Unmatched gen jets
271     t->Branch("nunmatch_j3",&jets_.nunmatch_j3,"nunmatch_j3/I");
272     t->Branch("unmatchpt_j3",jets_.unmatchpt_j3,"unmatchpt_j3[nunmatch_j3]/F");
273     t->Branch("unmatcheta_j3",jets_.unmatcheta_j3,"unmatcheta_j3[nunmatch_j3]/F");
274     t->Branch("unmatchy_j3",jets_.unmatchy_j3,"unmatchy_j3[nunmatch_j3]/F");
275     t->Branch("unmatchphi_j3",jets_.unmatchphi_j3,"unmatchphi_j3[nunmatch_j3]/F");
276     }
277    
278     // J4
279    
280     t->Branch("rawpt_j4",jets_.rawpt_j4,"rawpt_j4[nj4]/F");
281     t->Branch("jtpt_j4",jets_.jtpt_j4,"jtpt_j4[nj4]/F");
282     t->Branch("jteta_j4",jets_.jteta_j4,"jteta_j4[nj4]/F");
283     t->Branch("jty_j4",jets_.jty_j4,"jty_j4[nj4]/F");
284     t->Branch("jtphi_j4",jets_.jtphi_j4,"jtphi_j4[nj4]/F");
285     t->Branch("preL1et_j4",jets_.preL1et_j4,"preL1et_j4[nj4]/F");
286     t->Branch("L2_j4",jets_.L2_j4,"L2_j4[nj4]/F");
287     t->Branch("L3_j4",jets_.L3_j4,"L3_j4[nj4]/F");
288     t->Branch("area_j4",jets_.area_j4,"area_j4[nj4]/F");
289    
290     if(isMC_){
291     t->Branch("refpt_j4",jets_.refpt_j4,"refpt_j4[nj4]/F");
292     t->Branch("refeta_j4",jets_.refeta_j4,"refeta_j4[nj4]/F");
293     t->Branch("refy_j4",jets_.refy_j4,"refy_j4[nj4]/F");
294     t->Branch("refphi_j4",jets_.refphi_j4,"refphi_j4[nj4]/F");
295     t->Branch("refdrjt_j4",jets_.refdrjt_j4,"refdrjt_j4[nj4]/F");
296     t->Branch("refpartonpt_j4",jets_.refpartonpt_j4,"refpartonpt_j4[nj4]/F");
297     t->Branch("refpartonflavor_j4",jets_.refpartonflavor_j4,"refpartonflavor_j4[nj4]/F");
298    
299     // Unmatched gen jets
300     t->Branch("nunmatch_j4",&jets_.nunmatch_j4,"nunmatch_j4/I");
301     t->Branch("unmatchpt_j4",jets_.unmatchpt_j4,"unmatchpt_j4[nunmatch_j4]/F");
302     t->Branch("unmatcheta_j4",jets_.unmatcheta_j4,"unmatcheta_j4[nunmatch_j4]/F");
303     t->Branch("unmatchy_j4",jets_.unmatchy_j4,"unmatchy_j4[nunmatch_j4]/F");
304     t->Branch("unmatchphi_j4",jets_.unmatchphi_j4,"unmatchphi_j4[nunmatch_j4]/F");
305     }
306    
307     t->Branch("nPFcand",&jets_.nPFcand,"nPFcand/I");
308     t->Branch("candId",jets_.candId,"candId[nPFcand]/I");
309     t->Branch("candpt",jets_.candpt,"candpt[nPFcand]/F");
310     t->Branch("candeta",jets_.candeta,"candeta[nPFcand]/F");
311     //t->Branch("candy",jets_.candy,"candy[nPFcand]/F");
312     t->Branch("candphi",jets_.candphi,"candphi[nPFcand]/F");
313    
314    
315    
316     t->Branch("ntrack",&jets_.ntrack,"ntrack/I");
317     t->Branch("tracknhits",jets_.tracknhits,"tracknhits[ntrack]/I");
318     t->Branch("trackpt",jets_.trackpt,"trackpt[ntrack]/F");
319     t->Branch("tracketa",jets_.tracketa,"tracketa[ntrack]/F");
320     t->Branch("trackphi",jets_.trackphi,"trackphi[ntrack]/F");
321     t->Branch("tracksumecal",jets_.tracksumecal,"tracksumecal[ntrack]/F");
322     t->Branch("tracksumhcal",jets_.tracksumhcal,"tracksumhcal[ntrack]/F");
323     t->Branch("trackqual",jets_.trackqual,"trackqual[ntrack]/I");
324     t->Branch("chi2",jets_.trackchi2,"chi2[ntrack]/F");
325     t->Branch("chi2hit1D",jets_.trackchi2hit1D,"chi2hit1D[ntrack]/F");
326    
327     t->Branch("ptErr",jets_.trackptErr,"ptErr[ntrack]/F");
328    
329     t->Branch("d0Err",jets_.trackd0Err,"d0Err[ntrack]/F");
330     t->Branch("dzErr",jets_.trackdzErr,"dzErr[ntrack]/F");
331    
332     t->Branch("d0ErrTrk",jets_.trackd0ErrTrk,"d0ErrTrk[ntrack]/F");
333     t->Branch("dzErrTrk",jets_.trackdzErrTrk,"dzErrTrk[ntrack]/F");
334    
335     t->Branch("d0",jets_.trackd0,"d0[ntrack]/F");
336     t->Branch("dz",jets_.trackdz,"dz[ntrack]/F");
337    
338     // t->Branch("d0ErrBS",jets_.trackd0ErrBS,"d0ErrBS[ntrack]/F");
339     // t->Branch("dzErrBS",jets_.trackdzErrBS,"dzErrBS[ntrack]/F");
340     // t->Branch("d0BS",jets_.trackd0BS,"d0BS[ntrack]/F");
341     // t->Branch("dzBS",jets_.trackdzBS,"dzBS[ntrack]/F");
342    
343     t->Branch("nlayer",jets_.trackNlayer,"nlayer[ntrack]/I");
344     t->Branch("nlayer3D",jets_.trackNlayer3D,"nlayer3D[ntrack]/I");
345    
346     if(isMC_){
347     t->Branch("pthat",&jets_.pthat,"pthat/F");
348     t->Branch("trackfake",jets_.trackfake,"trackfake[ntrack]/I");
349     t->Branch("parton1_flavor",&jets_.parton1_flavor,"parton1_flavor/I");
350     t->Branch("parton2_flavor",&jets_.parton2_flavor,"parton2_flavor/I");
351     t->Branch("parton1_pt",&jets_.parton1_pt,"parton1_pt/F");
352     t->Branch("parton2_pt",&jets_.parton2_pt,"parton2_pt/F");
353     t->Branch("parton2_eta",&jets_.parton2_eta,"parton2_eta/F");
354     t->Branch("parton1_eta",&jets_.parton1_eta,"parton1_eta/F");
355     t->Branch("parton2_phi",&jets_.parton2_phi,"parton2_phi/F");
356     t->Branch("parton1_phi",&jets_.parton1_phi,"parton1_phi/F");
357     t->Branch("parton1_y",&jets_.parton1_y,"parton1_y/F");
358     t->Branch("parton2_y",&jets_.parton2_y,"parton2_y/F");
359     }
360     if(genParticleThresh_>0){
361     t->Branch("ngenp",&jets_.ngenp,"ngenp/I");
362     t->Branch("genppdgId",jets_.genppdgId,"genppdgId[ngenp]/I");
363     t->Branch("genppt",jets_.genppt,"genppt[ngenp]/F");
364     t->Branch("genpeta",jets_.genpeta,"genpeta[ngenp]/F");
365     t->Branch("genpphi",jets_.genpphi,"genpphi[ngenp]/F");
366     }
367     if(!isMC_){
368     t->Branch("nL1TBit",&jets_.nL1TBit,"nL1TBit/I");
369     t->Branch("l1TBit",jets_.l1TBit,"l1TBit[nL1TBit]/O");
370    
371     t->Branch("nL1ABit",&jets_.nL1ABit,"nL1ABit/I");
372     t->Branch("l1ABit",jets_.l1ABit,"l1ABit[nL1ABit]/O");
373    
374     t->Branch("nHLTBit",&jets_.nHLTBit,"nHLTBit/I");
375     t->Branch("hltBit",jets_.hltBit,"hltBit[nHLTBit]/O");
376    
377     }
378    
379    
380     TH1D::SetDefaultSumw2();
381    
382    
383     }
384    
385    
386     void
387     HiPFJetAnalyzer::analyze(const Event& iEvent,
388     const EventSetup& iSetup) {
389    
390    
391    
392     int event = iEvent.id().event();
393     int run = iEvent.id().run();
394     int lumi = iEvent.id().luminosityBlock();
395    
396     jets_.run = run;
397     jets_.evt = event;
398     jets_.lumi = lumi;
399    
400     LogDebug("HiPFJetAnalyzer")<<"START event: "<<event<<" in run "<<run<<endl;
401    
402     bool hasVertex = 0;
403     edm::Handle<reco::BeamSpot> beamSpotH;
404     iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpotH);
405    
406     edm::Handle<vector<reco::Vertex> >vertex;
407     iEvent.getByLabel(vertexTag_, vertex);
408    
409    
410    
411     if(vertex->size()>0 || vertex->begin()->isFake()) {
412     hasVertex = 1;
413     jets_.vx=vertex->begin()->x();
414     jets_.vy=vertex->begin()->y();
415     jets_.vz=vertex->begin()->z();
416    
417     jets_.vxErr = vertex->begin()->xError();
418     jets_.vyErr = vertex->begin()->yError();
419     jets_.vzErr = vertex->begin()->zError();
420     }else{
421     jets_.vx=beamSpotH->position().x();
422     jets_.vy=beamSpotH->position().y();
423     jets_.vz= 0;
424    
425     jets_.vxErr = beamSpotH->BeamWidthX();
426     jets_.vyErr = beamSpotH->BeamWidthY();
427     jets_.vzErr = 0;
428     }
429    
430     jets_.bsx=beamSpotH->position().x();
431     jets_.bsy=beamSpotH->position().y();
432     jets_.bsz= beamSpotH->position().z();
433     jets_.bswx=beamSpotH->BeamWidthX();
434     jets_.bswy=beamSpotH->BeamWidthY();
435    
436     jets_.hasVtx = hasVertex;
437    
438     int bin = -1;
439     double hf = 0.;
440     double b = 999.;
441    
442     if(useCentrality_){
443     //if(!isMC_){
444    
445     if(!centrality_) centrality_ = new CentralityProvider(iSetup);
446     centrality_->newEvent(iEvent,iSetup); // make sure you do this first in every event
447     //double c = centrality_->centralityValue();
448     const reco::Centrality *cent = centrality_->raw();
449    
450     hf = cent->EtHFhitSum();
451    
452     bin = centrality_->getBin();
453     b = centrality_->bMean();
454     //}
455     /*
456     else{
457     TFile * centFile = new TFile("/net/hidsk0001/d00/scratch/mnguyen/CMSSW_3_9_1_patch1/src/macros/Hydjet_CentTable.root");
458    
459     edm::Handle<reco::Centrality> cent;
460     iEvent.getByLabel(edm::InputTag("hiCentrality"),cent);
461     //cout<<" grabbed centrality "<<endl;
462     CentralityBins::RunMap cmap = getCentralityFromFile(centFile, "makeCentralityTableTFile", "HFhitsHydjet_2760GeV", 1, 1);
463    
464     // Still getting cent from hits. When tower tables become available, we need to switch
465     hf = cent->EtHFhitSum();
466     //cout<<" hf "<<hf<<endl;
467     bin = cmap[run]->getBin(hf);
468     b = cmap[run]->bMeanOfBin(bin);
469     }
470     */
471     }
472    
473    
474    
475     // not used, taking all jet
476     //double jetPtMin = 35.;
477    
478    
479     // loop the events
480    
481     jets_.bin = bin;
482     jets_.hf = hf;
483    
484     if(!isMC_){
485     fillL1Bits(iEvent);
486     fillHLTBits(iEvent);
487     }
488    
489     edm::Handle<pat::JetCollection> jets1;
490     iEvent.getByLabel(jetTag1_, jets1);
491    
492     edm::Handle<pat::JetCollection> jets2;
493     iEvent.getByLabel(jetTag2_, jets2);
494    
495     edm::Handle<pat::JetCollection> jets3;
496     iEvent.getByLabel(jetTag3_, jets3);
497    
498     edm::Handle<pat::JetCollection> jets4;
499     iEvent.getByLabel(jetTag4_, jets4);
500    
501    
502     edm::Handle< edm::View<reco::CaloJet> > recoJetColl1;
503     iEvent.getByLabel(recoJetTag1_, recoJetColl1 );
504    
505     edm::Handle< edm::View<reco::PFJet> > recoJetColl2;
506     iEvent.getByLabel(recoJetTag2_, recoJetColl2 );
507    
508     edm::Handle< edm::View<reco::PFJet> > recoJetColl3;
509     iEvent.getByLabel(recoJetTag3_, recoJetColl3 );
510    
511     edm::Handle< edm::View<reco::PFJet> > recoJetColl4;
512     iEvent.getByLabel(recoJetTag4_, recoJetColl4 );
513    
514    
515     Handle<PFCandidateCollection> pfCandidates;
516     iEvent.getByLabel(pfCandidatesTag_, pfCandidates);
517    
518    
519     Handle<vector<Track> > tracks;
520     iEvent.getByLabel(trackTag_, tracks);
521    
522     // do reco-to-sim association
523    
524     Handle<TrackingParticleCollection> TPCollectionHfake;
525     Handle<edm::View<reco::Track> > trackCollection;
526     ESHandle<TrackAssociatorBase> theAssociator;
527     const TrackAssociatorByHits *theAssociatorByHits;
528     reco::RecoToSimCollection recSimColl;
529    
530     if(hasSimInfo_) {
531     iEvent.getByLabel(simTracksTag_,TPCollectionHfake);
532     iEvent.getByLabel(trackTag_,trackCollection);
533     iSetup.get<TrackAssociatorRecord>().get("TrackAssociatorByHits",theAssociator);
534     theAssociatorByHits = (const TrackAssociatorByHits*) theAssociator.product();
535     recSimColl= theAssociatorByHits->associateRecoToSim(trackCollection,TPCollectionHfake,&iEvent);
536     }
537    
538    
539    
540     // FILL JRA TREE
541    
542     jets_.b = b;
543    
544    
545    
546     jets_.nj1 = 0;
547     jets_.nunmatch_j1 = 0;
548    
549    
550    
551     for(unsigned int j = 0 ; j < jets1->size(); ++j){
552     const pat::Jet& jet1 = (*jets1)[j];
553    
554     //cout<<" jet pt "<<jet1.pt()<<endl;
555     //if(jet1.pt() < jetPtMin) continue;
556     jets_.rawpt_j1[jets_.nj1]=jet1.correctedJet("Uncorrected").pt();
557     jets_.jtpt_j1[jets_.nj1] = jet1.pt();
558     jets_.jteta_j1[jets_.nj1] = jet1.eta();
559     jets_.jtphi_j1[jets_.nj1] = jet1.phi();
560     jets_.jty_j1[jets_.nj1] = jet1.eta();
561    
562    
563     // cout<<" abs corr "<<jet1.corrFactor("L3Absolute")<<endl;
564     //cout<<" abs corr "<<jet1.corrFactor("L3Absolute")<<endl;
565    
566    
567     float L2Corr = jet1.correctedJet("L2Relative").pt()/jet1.correctedJet("Uncorrected").pt();
568     float L3Corr = jet1.correctedJet("L3Absolute").pt()/jet1.correctedJet("L2Relative").pt();
569    
570    
571     jets_.L2_j1[jets_.nj1] = L2Corr;
572     jets_.L3_j1[jets_.nj1] = L3Corr;
573    
574    
575     jets_.area_j1[jets_.nj1] = jet1.jetArea();
576    
577     // Match to reco jet to find unsubtracted jet energy
578    
579     if(1==0){
580     int recoJetSize = recoJetColl1->size();
581    
582     jets_.preL1et_j1[jets_.nj1] = -1;
583    
584     //cout<<" patJet_eta "<<jet1.eta()<<" patJet_phi "<<jet1.phi()<<" patJet_et "<<jet1.et()<<endl;
585    
586     for(int iRecoJet = 0; iRecoJet < recoJetSize; ++iRecoJet){
587    
588     reco::CaloJet recoJet1 = ((*recoJetColl1)[iRecoJet]);
589    
590    
591     double recoJet_eta = recoJet1.eta();
592     double recoJet_phi = recoJet1.phi();
593     //cout<<" recoJet_eta "<<recoJet_eta<<" recoJet_phi "<<recoJet_phi<<" recoJet_et "<<recoJet1.et()<<endl;
594    
595    
596     if(fabs(recoJet1.eta()-jet1.eta()) < 0.001
597     && fabs(acos(cos((recoJet1.phi()-jet1.phi())))) < 0.001)
598     {
599     jets_.preL1et_j1[jets_.nj1] = recoJet1.et();
600    
601     //cout<<"Match found, recoJet1.et "<<recoJet1.et()<< " recoJet1.eta "<<jet1.eta()<<" recoJet1.phi "<<recoJet1.phi()<<endl;
602     break;
603     }
604     }
605     if(jets_.preL1et_j1[jets_.nj1] == -1){
606    
607    
608     //There's a known issue here. If the background subtraction oversubtracts I've set the patJet.et() to zero. That would be fine if I had also set the eta and phi. We could then recover the pre-subtracted energy. However, I never bothered to set the eta and phi for theses jets (doh!). Next time I repass the data I won't be so stupid.
609    
610    
611    
612     if(jet1.et()>0)cout<<"Match *NOT* found, patJet1.et "<<jet1.et()<< " patJet1.eta "<<jet1.eta()<<" patJet1.phi() "<<jet1.phi()<<endl;
613     }
614     }
615     if(isMC_){
616    
617    
618     if(jet1.genJet()!=0 && jet1.genJet()->pt()>1.0 && jet1.genJet()->pt()<999999){
619     jets_.refpt_j1[jets_.nj1] = jet1.genJet()->pt();
620     jets_.refeta_j1[jets_.nj1] = jet1.genJet()->eta();
621     jets_.refphi_j1[jets_.nj1] = jet1.genJet()->phi();
622     jets_.refy_j1[jets_.nj1] = jet1.genJet()->eta();
623    
624     jets_.refdrjt_j1[jets_.nj1] = reco::deltaR(jet1,*(jet1.genJet()));
625     }
626     else{
627     jets_.refpt_j1[jets_.nj1] = 0;
628     jets_.refeta_j1[jets_.nj1] = -999;
629     jets_.refphi_j1[jets_.nj1] = -999;
630     jets_.refy_j1[jets_.nj1] = -999;
631     }
632    
633     if (jet1.genParton()) {
634     jets_.refpartonpt_j1[jets_.nj1] = jet1.genParton()->pt();
635     jets_.refpartonflavor_j1[jets_.nj1] = jet1.genParton()->pdgId();
636     } else {
637     jets_.refpartonpt_j1[jets_.nj1] = -999;
638     jets_.refpartonflavor_j1[jets_.nj1] = -999;
639     }
640    
641     }
642     jets_.nj1++;
643     }
644    
645     if(isMC_){
646     edm::Handle<vector<reco::GenJet> >genjets1;
647     iEvent.getByLabel(genJetTag1_, genjets1);
648    
649     for(unsigned int igen = 0 ; igen < genjets1->size(); ++igen){
650     const reco::GenJet & genjet1 = (*genjets1)[igen];
651    
652     float genjet_pt = genjet1.pt();
653    
654     // threshold to reduce size of output in minbias PbPb
655     if(genjet_pt>20.){
656    
657     int isMatched=0;
658    
659     for(unsigned int ijet = 0 ; ijet < jets1->size(); ++ijet){
660     const pat::Jet& jet1 = (*jets1)[ijet];
661    
662     if(jet1.genJet()){
663     if(fabs(genjet1.pt()-jet1.genJet()->pt())<0.0001 &&
664     fabs(genjet1.eta()-jet1.genJet()->eta())<0.0001 &&
665     (fabs(genjet1.phi()-jet1.genJet()->phi())<0.0001 || fabs(fabs(genjet1.phi()-jet1.genJet()->phi()) - 2.0*TMath::Pi()) < 0.0001 )){
666    
667     isMatched =1;
668     break;
669     }
670     }
671     }
672    
673     if(!isMatched){
674     jets_.unmatchpt_j1[jets_.nunmatch_j1] = genjet_pt;
675     jets_.unmatcheta_j1[jets_.nunmatch_j1] = genjet1.eta();
676     jets_.unmatchphi_j1[jets_.nunmatch_j1] = genjet1.phi();
677     jets_.unmatchy_j1[jets_.nunmatch_j1] = genjet1.eta();
678    
679     jets_.nunmatch_j1++;
680    
681     }
682    
683     }
684     }
685     }
686    
687    
688    
689    
690    
691     jets_.nj2 = 0;
692     jets_.nunmatch_j2 = 0;
693    
694    
695     for(unsigned int j = 0 ; j < jets2->size(); ++j){
696     const pat::Jet& jet2 = (*jets2)[j];
697    
698     //cout<<" jet pt "<<jet2.pt()<<endl;
699     //if(jet2.pt() < jetPtMin) continue;
700     jets_.rawpt_j2[jets_.nj2]=jet2.correctedJet("Uncorrected").pt();
701     jets_.jtpt_j2[jets_.nj2] = jet2.pt();
702     jets_.jteta_j2[jets_.nj2] = jet2.eta();
703     jets_.jtphi_j2[jets_.nj2] = jet2.phi();
704     jets_.jty_j2[jets_.nj2] = jet2.eta();
705     // cout<<" abs corr "<<jet2.corrFactor("L3Absolute")<<endl;
706     //cout<<" abs corr "<<jet2.corrFactor("L3Absolute")<<endl;
707    
708    
709     float L2Corr = jet2.correctedJet("L2Relative").pt()/jet2.correctedJet("Uncorrected").pt();
710     float L3Corr = jet2.correctedJet("L3Absolute").pt()/jet2.correctedJet("L2Relative").pt();
711    
712    
713     jets_.L2_j2[jets_.nj2] = L2Corr;
714     jets_.L3_j2[jets_.nj2] = L3Corr;
715    
716     jets_.area_j2[jets_.nj2] = jet2.jetArea();
717    
718     // Match to reco jet to find unsubtracted jet energy
719     if(1==0){
720     int recoJetSize2 = recoJetColl2->size();
721    
722     jets_.preL1et_j2[jets_.nj2] = -1;
723    
724     //cout<<" patJet_eta "<<jet2.eta()<<" patJet_phi "<<jet2.phi()<<" patJet_et "<<jet2.et()<<endl;
725    
726     for(int iRecoJet = 0; iRecoJet < recoJetSize2; ++iRecoJet){
727    
728     reco::PFJet recoJet2 = ((*recoJetColl2)[iRecoJet]);
729    
730    
731     double recoJet_eta = recoJet2.eta();
732     double recoJet_phi = recoJet2.phi();
733     //cout<<" recoJet_eta "<<recoJet_eta<<" recoJet_phi "<<recoJet_phi<<" recoJet_et "<<recoJet2.et()<<endl;
734    
735    
736     if(fabs(recoJet2.eta()-jet2.eta()) < 0.001
737     && fabs(acos(cos((recoJet2.phi()-jet2.phi())))) < 0.001)
738     {
739     jets_.preL1et_j2[jets_.nj2] = recoJet2.et();
740    
741     //cout<<"Match found, recoJet2.et "<<recoJet2.et()<< " recoJet2.eta "<<jet2.eta()<<" recoJet2.phi "<<recoJet2.phi()<<endl;
742     break;
743     }
744     }
745     if(jets_.preL1et_j2[jets_.nj2] == -1){
746    
747    
748     //There's a known issue here. If the background subtraction oversubtracts I've set the patJet.et() to zero. That would be fine if I had also set the eta and phi. We could then recover the pre-subtracted energy. However, I never bothered to set the eta and phi for theses jets (doh!). Next time I repass the data I won't be so stupid.
749    
750    
751    
752     if(jet2.et()>0)cout<<"Match *NOT* found, patJet2.et "<<jet2.et()<< " patJet2.eta "<<jet2.eta()<<" patJet2.phi() "<<jet2.phi()<<endl;
753     }
754     }
755     if(isMC_){
756    
757    
758     if(jet2.genJet()!=0 && jet2.genJet()->pt()>1.0 && jet2.genJet()->pt()<999999){
759     jets_.refpt_j2[jets_.nj2] = jet2.genJet()->pt();
760     jets_.refeta_j2[jets_.nj2] = jet2.genJet()->eta();
761     jets_.refphi_j2[jets_.nj2] = jet2.genJet()->phi();
762     jets_.refy_j2[jets_.nj2] = jet2.genJet()->eta();
763    
764     jets_.refdrjt_j2[jets_.nj2] = reco::deltaR(jet2,*(jet2.genJet()));
765     }
766     else{
767     jets_.refpt_j2[jets_.nj2] = 0;
768     jets_.refeta_j2[jets_.nj2] = -999;
769     jets_.refphi_j2[jets_.nj2] = -999;
770     jets_.refy_j2[jets_.nj2] = -999;
771     }
772    
773     if (jet2.genParton()) {
774     jets_.refpartonpt_j2[jets_.nj2] = jet2.genParton()->pt();
775     jets_.refpartonflavor_j2[jets_.nj2] = jet2.genParton()->pdgId();
776     } else {
777     jets_.refpartonpt_j2[jets_.nj2] = -999;
778     jets_.refpartonflavor_j2[jets_.nj2] = -999;
779     }
780     }
781    
782    
783     jets_.nj2++;
784    
785     }
786    
787     if(isMC_){
788    
789     edm::Handle<vector<reco::GenJet> >genjets2;
790     iEvent.getByLabel(genJetTag2_, genjets2);
791    
792     for(unsigned int igen = 0 ; igen < genjets2->size(); ++igen){
793     const reco::GenJet & genjet2 = (*genjets2)[igen];
794    
795     float genjet_pt = genjet2.pt();
796    
797     // threshold to reduce size of output in minbias PbPb
798     if(genjet_pt>20.){
799    
800     int isMatched=0;
801    
802     for(unsigned int ijet = 0 ; ijet < jets2->size(); ++ijet){
803     const pat::Jet& jet2 = (*jets2)[ijet];
804    
805     if(jet2.genJet()){
806     if(fabs(genjet2.pt()-jet2.genJet()->pt())<0.0001 &&
807     fabs(genjet2.eta()-jet2.genJet()->eta())<0.0001 &&
808     (fabs(genjet2.phi()-jet2.genJet()->phi())<0.0001 || fabs(fabs(genjet2.phi()-jet2.genJet()->phi()) - 2.0*TMath::Pi()) < 0.0001 )){
809    
810     isMatched =1;
811     break;
812     }
813     }
814     }
815    
816     if(!isMatched){
817     jets_.unmatchpt_j2[jets_.nunmatch_j2] = genjet_pt;
818     jets_.unmatcheta_j2[jets_.nunmatch_j2] = genjet2.eta();
819     jets_.unmatchphi_j2[jets_.nunmatch_j2] = genjet2.phi();
820     jets_.unmatchy_j2[jets_.nunmatch_j2] = genjet2.eta();
821    
822     jets_.nunmatch_j2++;
823    
824     }
825    
826     }
827     }
828     }
829    
830    
831    
832     jets_.nj3 = 0;
833     jets_.nunmatch_j3 = 0;
834    
835     //cout<<" jets size "<<jets->size()<<endl;
836    
837    
838     for(unsigned int j = 0 ; j < jets3->size(); ++j){
839     const pat::Jet& jet3 = (*jets3)[j];
840    
841     //cout<<" jet pt "<<jet3.pt()<<endl;
842     //if(jet3.pt() < jetPtMin) continue;
843     jets_.rawpt_j3[jets_.nj3]=jet3.correctedJet("Uncorrected").pt();
844     jets_.jtpt_j3[jets_.nj3] = jet3.pt();
845     jets_.jteta_j3[jets_.nj3] = jet3.eta();
846     jets_.jtphi_j3[jets_.nj3] = jet3.phi();
847     jets_.jty_j3[jets_.nj3] = jet3.eta();
848     // cout<<" abs corr "<<jet3.corrFactor("L3Absolute")<<endl;
849     //cout<<" abs corr "<<jet3.corrFactor("L3Absolute")<<endl;
850    
851    
852    
853     float L2Corr = jet3.correctedJet("L2Relative").pt()/jet3.correctedJet("Uncorrected").pt();
854     float L3Corr = jet3.correctedJet("L3Absolute").pt()/jet3.correctedJet("L2Relative").pt();
855    
856    
857     jets_.L2_j3[jets_.nj3] = L2Corr;
858     jets_.L3_j3[jets_.nj3] = L3Corr;
859    
860     jets_.area_j3[jets_.nj3] = jet3.jetArea();
861    
862     // Match to reco jet to find unsubtracted jet energy
863     if(1==0){
864     int recoJetSize3 = recoJetColl3->size();
865    
866     jets_.preL1et_j3[jets_.nj3] = -1;
867    
868     //cout<<" patJet_eta "<<jet3.eta()<<" patJet_phi "<<jet3.phi()<<" patJet_et "<<jet3.et()<<endl;
869    
870     for(int iRecoJet = 0; iRecoJet < recoJetSize3; ++iRecoJet){
871    
872     reco::PFJet recoJet3 = ((*recoJetColl3)[iRecoJet]);
873    
874    
875     double recoJet_eta = recoJet3.eta();
876     double recoJet_phi = recoJet3.phi();
877     //cout<<" recoJet_eta "<<recoJet_eta<<" recoJet_phi "<<recoJet_phi<<" recoJet_et "<<recoJet3.et()<<endl;
878    
879    
880     if(fabs(recoJet3.eta()-jet3.eta()) < 0.001
881     && fabs(acos(cos((recoJet3.phi()-jet3.phi())))) < 0.001)
882     {
883     jets_.preL1et_j3[jets_.nj3] = recoJet3.et();
884    
885     //cout<<"Match found, recoJet3.et "<<recoJet3.et()<< " recoJet3.eta "<<jet3.eta()<<" recoJet3.phi "<<recoJet3.phi()<<endl;
886     break;
887     }
888     }
889     if(jets_.preL1et_j3[jets_.nj3] == -1){
890    
891    
892     // There's a known issue here. If the background subtraction oversubtracts I've set the patJet.et() to zero. That would be fine if I had also set the eta and phi. We could then recover the pre-subtracted energy. However, I never bothered to set the eta and phi for theses jets (doh!). Next time I repass the data I won't be so stupid.
893    
894    
895    
896     if(jet3.et()>0)cout<<"Match *NOT* found, patJet3.et "<<jet3.et()<< " patJet3.eta "<<jet3.eta()<<" patJet3.phi() "<<jet3.phi()<<endl;
897     }
898     }
899     if(isMC_){
900    
901    
902     if(jet3.genJet()!=0 && jet3.genJet()->pt()>1.0 && jet3.genJet()->pt()<999999){
903     jets_.refpt_j3[jets_.nj3] = jet3.genJet()->pt();
904     jets_.refeta_j3[jets_.nj3] = jet3.genJet()->eta();
905     jets_.refphi_j3[jets_.nj3] = jet3.genJet()->phi();
906     jets_.refy_j3[jets_.nj3] = jet3.genJet()->eta();
907    
908     jets_.refdrjt_j3[jets_.nj3] = reco::deltaR(jet3,*(jet3.genJet()));
909     }
910     else{
911     jets_.refpt_j3[jets_.nj3] = 0;
912     jets_.refeta_j3[jets_.nj3] = -999;
913     jets_.refphi_j3[jets_.nj3] = -999;
914     jets_.refy_j3[jets_.nj3] = -999;
915     }
916    
917     if (jet3.genParton()) {
918     jets_.refpartonpt_j3[jets_.nj3] = jet3.genParton()->pt();
919     jets_.refpartonflavor_j3[jets_.nj3] = jet3.genParton()->pdgId();
920     } else {
921     jets_.refpartonpt_j3[jets_.nj3] = -999;
922     jets_.refpartonflavor_j3[jets_.nj3] = -999;
923     }
924    
925     }
926    
927    
928    
929     jets_.nj3++;
930    
931    
932     }
933    
934     if(isMC_){
935    
936     edm::Handle<vector<reco::GenJet> >genjets3;
937     iEvent.getByLabel(genJetTag3_, genjets3);
938    
939     for(unsigned int igen = 0 ; igen < genjets3->size(); ++igen){
940     const reco::GenJet & genjet3 = (*genjets3)[igen];
941    
942     float genjet_pt = genjet3.pt();
943    
944     // threshold to reduce size of output in minbias PbPb
945     if(genjet_pt>20.){
946    
947     int isMatched=0;
948    
949     for(unsigned int ijet = 0 ; ijet < jets3->size(); ++ijet){
950     const pat::Jet& jet3 = (*jets3)[ijet];
951    
952     if(jet3.genJet()){
953     if(fabs(genjet3.pt()-jet3.genJet()->pt())<0.0001 &&
954     fabs(genjet3.eta()-jet3.genJet()->eta())<0.0001 &&
955     (fabs(genjet3.phi()-jet3.genJet()->phi())<0.0001 || fabs(fabs(genjet3.phi()-jet3.genJet()->phi()) - 2.0*TMath::Pi()) < 0.0001 )){
956    
957     isMatched =1;
958     break;
959     }
960     }
961     }
962    
963     if(!isMatched){
964     jets_.unmatchpt_j3[jets_.nunmatch_j3] = genjet_pt;
965     jets_.unmatcheta_j3[jets_.nunmatch_j3] = genjet3.eta();
966     jets_.unmatchphi_j3[jets_.nunmatch_j3] = genjet3.phi();
967     jets_.unmatchy_j3[jets_.nunmatch_j3] = genjet3.eta();
968    
969     jets_.nunmatch_j3++;
970    
971     }
972    
973     }
974     }
975     }
976    
977    
978    
979    
980     jets_.nj4 = 0;
981     jets_.nunmatch_j4 = 0;
982    
983    
984     for(unsigned int j = 0 ; j < jets4->size(); ++j){
985     const pat::Jet& jet4 = (*jets4)[j];
986    
987     //cout<<" jet pt "<<jet4.pt()<<endl;
988     //if(jet4.pt() < jetPtMin) continue;
989     jets_.rawpt_j4[jets_.nj4]=jet4.correctedJet("Uncorrected").pt();
990     jets_.jtpt_j4[jets_.nj4] = jet4.pt();
991     jets_.jteta_j4[jets_.nj4] = jet4.eta();
992     jets_.jtphi_j4[jets_.nj4] = jet4.phi();
993     jets_.jty_j4[jets_.nj4] = jet4.eta();
994     // cout<<" abs corr "<<jet4.corrFactor("L3Absolute")<<endl;
995     //cout<<" abs corr "<<jet4.corrFactor("L3Absolute")<<endl;
996    
997    
998     float L2Corr = jet4.correctedJet("L2Relative").pt()/jet4.correctedJet("Uncorrected").pt();
999     float L3Corr = jet4.correctedJet("L3Absolute").pt()/jet4.correctedJet("L2Relative").pt();
1000    
1001    
1002     jets_.L2_j4[jets_.nj4] = L2Corr;
1003     jets_.L3_j4[jets_.nj4] = L3Corr;
1004    
1005     jets_.area_j4[jets_.nj4] = jet4.jetArea();
1006    
1007     // Match to reco jet to find unsubtracted jet energy
1008     if(1==0){
1009     int recoJetSize4 = recoJetColl4->size();
1010    
1011     jets_.preL1et_j4[jets_.nj4] = -1;
1012    
1013     //cout<<" patJet_eta "<<jet4.eta()<<" patJet_phi "<<jet4.phi()<<" patJet_et "<<jet4.et()<<endl;
1014    
1015     for(int iRecoJet = 0; iRecoJet < recoJetSize4; ++iRecoJet){
1016    
1017     reco::PFJet recoJet4 = ((*recoJetColl4)[iRecoJet]);
1018    
1019    
1020     double recoJet_eta = recoJet4.eta();
1021     double recoJet_phi = recoJet4.phi();
1022     //cout<<" recoJet_eta "<<recoJet_eta<<" recoJet_phi "<<recoJet_phi<<" recoJet_et "<<recoJet4.et()<<endl;
1023    
1024    
1025     if(fabs(recoJet4.eta()-jet4.eta()) < 0.001
1026     && fabs(acos(cos((recoJet4.phi()-jet4.phi())))) < 0.001)
1027     {
1028     jets_.preL1et_j4[jets_.nj4] = recoJet4.et();
1029    
1030     //cout<<"Match found, recoJet4.et "<<recoJet4.et()<< " recoJet4.eta "<<jet4.eta()<<" recoJet4.phi "<<recoJet4.phi()<<endl;
1031     break;
1032     }
1033     }
1034     if(jets_.preL1et_j4[jets_.nj4] == -1){
1035    
1036    
1037     //There's a known issue here. If the background subtraction oversubtracts I've set the patJet.et() to zero. That would be fine if I had also set the eta and phi. We could then recover the pre-subtracted energy. However, I never bothered to set the eta and phi for theses jets (doh!). Next time I repass the data I won't be so stupid.
1038    
1039    
1040    
1041     if(jet4.et()>0)cout<<"Match *NOT* found, patJet4.et "<<jet4.et()<< " patJet4.eta "<<jet4.eta()<<" patJet4.phi() "<<jet4.phi()<<endl;
1042     }
1043     }
1044     if(isMC_){
1045    
1046    
1047     if(jet4.genJet()!=0 && jet4.genJet()->pt()>1.0 && jet4.genJet()->pt()<999999){
1048     jets_.refpt_j4[jets_.nj4] = jet4.genJet()->pt();
1049     jets_.refeta_j4[jets_.nj4] = jet4.genJet()->eta();
1050     jets_.refphi_j4[jets_.nj4] = jet4.genJet()->phi();
1051     jets_.refy_j4[jets_.nj4] = jet4.genJet()->eta();
1052    
1053     jets_.refdrjt_j4[jets_.nj4] = reco::deltaR(jet4,*(jet4.genJet()));
1054     }
1055     else{
1056     jets_.refpt_j4[jets_.nj4] = 0;
1057     jets_.refeta_j4[jets_.nj4] = -999;
1058     jets_.refphi_j4[jets_.nj4] = -999;
1059     jets_.refy_j4[jets_.nj4] = -999;
1060     }
1061    
1062     if (jet4.genParton()) {
1063     jets_.refpartonpt_j4[jets_.nj4] = jet4.genParton()->pt();
1064     jets_.refpartonflavor_j4[jets_.nj4] = jet4.genParton()->pdgId();
1065     } else {
1066     jets_.refpartonpt_j4[jets_.nj4] = -999;
1067     jets_.refpartonflavor_j4[jets_.nj4] = -999;
1068     }
1069    
1070     }
1071    
1072    
1073     jets_.nj4++;
1074    
1075     }
1076    
1077     if(isMC_){
1078    
1079     edm::Handle<vector<reco::GenJet> >genjets4;
1080     iEvent.getByLabel(genJetTag4_, genjets4);
1081    
1082     for(unsigned int igen = 0 ; igen < genjets4->size(); ++igen){
1083     const reco::GenJet & genjet4 = (*genjets4)[igen];
1084    
1085     float genjet_pt = genjet4.pt();
1086    
1087     // threshold to reduce size of output in minbias PbPb
1088     if(genjet_pt>20.){
1089    
1090     int isMatched=0;
1091    
1092     for(unsigned int ijet = 0 ; ijet < jets4->size(); ++ijet){
1093     const pat::Jet& jet4 = (*jets4)[ijet];
1094    
1095     if(jet4.genJet()){
1096     if(fabs(genjet4.pt()-jet4.genJet()->pt())<0.0001 &&
1097     fabs(genjet4.eta()-jet4.genJet()->eta())<0.0001 &&
1098     (fabs(genjet4.phi()-jet4.genJet()->phi())<0.0001 || fabs(fabs(genjet4.phi()-jet4.genJet()->phi()) - 2.0*TMath::Pi()) < 0.0001 )){
1099    
1100     isMatched =1;
1101     break;
1102     }
1103     }
1104     }
1105    
1106     if(!isMatched){
1107     jets_.unmatchpt_j4[jets_.nunmatch_j4] = genjet_pt;
1108     jets_.unmatcheta_j4[jets_.nunmatch_j4] = genjet4.eta();
1109     jets_.unmatchphi_j4[jets_.nunmatch_j4] = genjet4.phi();
1110     jets_.unmatchy_j4[jets_.nunmatch_j4] = genjet4.eta();
1111    
1112     jets_.nunmatch_j4++;
1113    
1114     }
1115    
1116     }
1117     }
1118     }
1119    
1120    
1121    
1122     for( unsigned icand=0; icand<pfCandidates->size(); icand++ ) {
1123    
1124     const reco::PFCandidate& cand = (*pfCandidates)[icand];
1125    
1126     float particleEta = cand.eta();
1127    
1128     //if(fabs(particleEta)>2.5) continue;
1129    
1130    
1131    
1132     // PF PId Convention:
1133     // 1 = Charged Hadrons
1134     // 2 = Electrons (not included)
1135     // 3 = Muons
1136     // 4 = Photons
1137     // 5 = Neutral Hadrons
1138    
1139    
1140    
1141    
1142     int particleId = (int)cand.particleId();
1143     float particlePt = cand.pt();
1144    
1145     if(particlePt<0.3) continue;
1146    
1147    
1148     // can use varid thresholds if we want
1149     //if(particleId==1 && particlePt < 0.9) continue;
1150     //if(particleId==3 && particlePt < 0.9) continue;
1151     //if(particleId==4 && particlePt < 0.3) continue;
1152     //if(particleId==5 && particlePt < 0.9) continue;
1153    
1154    
1155     if(particleId==3&&particlePt>100) cout<<" likely a badly reconstructed MUON "<<endl;
1156    
1157     jets_.candId[jets_.nPFcand] = particleId;
1158     jets_.candpt[jets_.nPFcand] = particlePt;
1159     jets_.candeta[jets_.nPFcand] = particleEta;
1160     jets_.candphi[jets_.nPFcand] = cand.phi();
1161     //jets_.candy[jets_.nPFcand] = cand.y();
1162    
1163     if(particleId==3&&particlePt>100) cout<<" found a misreconstructed MUON, pT = "<<particlePt<<endl;
1164    
1165     jets_.nPFcand++;
1166    
1167     //cout<<" jets_.nPFcand "<<jets_.nPFcand<<endl;
1168     }
1169    
1170     //cout<<" ntracks: "<<tracks->size()<<endl;
1171    
1172     for(unsigned int it=0; it<tracks->size(); ++it){
1173     const reco::Track & track = (*tracks)[it];
1174    
1175     int count1dhits = 0;
1176     double chi2n_hit1D = 0;
1177     trackingRecHit_iterator edh = track.recHitsEnd();
1178     for (trackingRecHit_iterator ith = track.recHitsBegin(); ith != edh; ++ith) {
1179     const TrackingRecHit * hit = ith->get();
1180     DetId detid = hit->geographicalId();
1181     if (hit->isValid()) {
1182     if (typeid(*hit) == typeid(SiStripRecHit1D)) ++count1dhits;
1183     }
1184     }
1185     if (count1dhits > 0) {
1186     double chi2 = track.chi2();
1187     double ndof = track.ndof();
1188     chi2n_hit1D = (chi2+count1dhits)/double(ndof+count1dhits);
1189     }
1190    
1191     // Could makes some track selection here
1192     jets_.tracknhits[jets_.ntrack] = track.numberOfValidHits();
1193     jets_.trackpt[jets_.ntrack] = track.pt();
1194     jets_.tracketa[jets_.ntrack] = track.eta();
1195     jets_.trackphi[jets_.ntrack] = track.phi();
1196    
1197     jets_.trackptErr[jets_.ntrack] = track.ptError();
1198     jets_.trackchi2[jets_.ntrack] = track.normalizedChi2();
1199     jets_.trackchi2hit1D[jets_.ntrack] = chi2n_hit1D;
1200    
1201     jets_.tracksumecal[jets_.ntrack] = 0.;
1202     jets_.tracksumhcal[jets_.ntrack] = 0.;
1203    
1204     reco::TrackBase::TrackQuality trackQualityTight = TrackBase::qualityByName("highPurity");
1205     jets_.trackqual[jets_.ntrack]=(int)track.quality(trackQualityTight);
1206    
1207     jets_.trackfake[jets_.ntrack]=0;
1208    
1209     reco::TrackRef trackRef=reco::TrackRef(tracks,it);
1210    
1211     if(hasVertex){
1212     jets_.trackd0[jets_.ntrack] = -track.dxy(vertex->begin()->position());
1213     jets_.trackdz[jets_.ntrack] = track.dz(vertex->begin()->position());
1214     jets_.trackd0Err[jets_.ntrack] = sqrt ( (track.d0Error()*track.d0Error()) + (vertex->begin()->xError()*vertex->begin()->yError()) );
1215     jets_.trackdzErr[jets_.ntrack] = sqrt ( (track.dzError()*track.dzError()) + (vertex->begin()->zError()*vertex->begin()->zError()) );
1216     }else{
1217     jets_.trackd0[jets_.ntrack] = -track.dxy(beamSpotH->position());
1218     jets_.trackdz[jets_.ntrack] = 0;
1219     jets_.trackd0Err[jets_.ntrack] = sqrt ( (track.d0Error()*track.d0Error()) + (beamSpotH->BeamWidthX()*beamSpotH->BeamWidthY()) );
1220     jets_.trackdzErr[jets_.ntrack] = 0;
1221     }
1222    
1223     jets_.trackd0BS[jets_.ntrack] = -track.dxy(beamSpotH->position());
1224     jets_.trackdzBS[jets_.ntrack] = track.dz(beamSpotH->position());
1225     jets_.trackd0ErrBS[jets_.ntrack] = sqrt ( (track.d0Error()*track.d0Error()) + (beamSpotH->BeamWidthX()*beamSpotH->BeamWidthY()) );
1226     jets_.trackdzErrBS[jets_.ntrack] = 0;
1227    
1228     jets_.trackd0ErrTrk[jets_.ntrack] = track.d0Error();
1229     jets_.trackdzErrTrk[jets_.ntrack] = track.dzError();
1230    
1231     jets_.trackNlayer[jets_.ntrack] = track.hitPattern().trackerLayersWithMeasurement();
1232     jets_.trackNlayer3D[jets_.ntrack] = track.hitPattern().pixelLayersWithMeasurement() + track.hitPattern().numberOfValidStripLayersWithMonoAndStereo();
1233    
1234     if(hasSimInfo_)
1235     if(recSimColl.find(edm::RefToBase<reco::Track>(trackRef)) == recSimColl.end())
1236     jets_.trackfake[jets_.ntrack]=1;
1237    
1238    
1239    
1240     int pfCandMatchFound = 0;
1241    
1242     // loop over pf candidates to get calo-track matching info
1243     for( unsigned icand=0; icand<pfCandidates->size(); icand++ ) {
1244    
1245     const reco::PFCandidate& cand = (*pfCandidates)[icand];
1246    
1247     float cand_type = cand.particleId();
1248    
1249     // only charged hadrons and leptons can be asscociated with a track
1250     if(!(cand_type == PFCandidate::h || //type1
1251     cand_type == PFCandidate::e || //type2
1252     cand_type == PFCandidate::mu //type3
1253     )
1254     ) continue;
1255    
1256    
1257     // if working with 2 different track collections this doesn't work
1258     if(cand.trackRef() != trackRef) continue;
1259     //if(fabs(cand.pt()-track.pt())>0.001||fabs(cand.eta()-track.eta())>0.001||fabs(acos(cos(cand.phi()-track.phi())))>0.001) continue;
1260    
1261     pfCandMatchFound = 1;
1262    
1263     for(unsigned iblock=0; iblock<cand.elementsInBlocks().size(); iblock++) {
1264    
1265     PFBlockRef blockRef = cand.elementsInBlocks()[iblock].first;
1266     unsigned indexInBlock = cand.elementsInBlocks()[iblock].second;
1267    
1268    
1269     const edm::OwnVector< reco::PFBlockElement>& elements = (*blockRef).elements();
1270    
1271     //This tells you what type of element it is:
1272     //cout<<" block type"<<elements[indexInBlock].type()<<endl;
1273    
1274     switch (elements[indexInBlock].type()) {
1275    
1276     case PFBlockElement::ECAL: {
1277     reco::PFClusterRef clusterRef = elements[indexInBlock].clusterRef();
1278     double eet = clusterRef->energy()/cosh(clusterRef->eta());
1279     if(verbose_)cout<<" ecal energy "<<clusterRef->energy()<<endl;
1280     jets_.tracksumecal[jets_.ntrack] += eet;
1281     break;
1282     }
1283    
1284     case PFBlockElement::HCAL: {
1285     reco::PFClusterRef clusterRef = elements[indexInBlock].clusterRef();
1286     double eet = clusterRef->energy()/cosh(clusterRef->eta());
1287     if(verbose_)cout<<" hcal energy "<<clusterRef->energy()<<endl;
1288     jets_.tracksumhcal[jets_.ntrack] += eet;
1289     break;
1290     }
1291     case PFBlockElement::TRACK: {
1292     //This is just the reference to the track itself, since tracks can never be linked to other tracks
1293     break;
1294     }
1295     default:
1296     break;
1297     }
1298     // Could do more stuff here, e.g., pre-shower, HF
1299    
1300     }
1301    
1302     }
1303    
1304     if(!pfCandMatchFound){
1305     jets_.tracksumecal[jets_.ntrack] =-1;
1306     jets_.tracksumhcal[jets_.ntrack] =-1;
1307    
1308     }
1309    
1310     jets_.ntrack++;
1311    
1312     }
1313    
1314     // make configurable, so that gen particles aren't run with MB
1315    
1316     if(isMC_){
1317    
1318     edm::Handle<GenEventInfoProduct> hEventInfo;
1319     iEvent.getByLabel(eventInfoTag_,hEventInfo);
1320    
1321     jets_.pthat = hEventInfo->qScale();
1322    
1323     //getPartons(iEvent, iSetup );
1324    
1325     if(genParticleThresh_>0){
1326     edm::Handle <reco::GenParticleCollection> genParticles;
1327     iEvent.getByLabel (genParticleTag_, genParticles );
1328    
1329    
1330     for( unsigned igen=0; igen<genParticles->size(); igen++ ) {
1331    
1332    
1333     const reco::GenParticle & genp = (*genParticles)[igen];
1334    
1335     if(genp.status()!=1) continue;
1336    
1337     jets_.genppt[jets_.ngenp] = genp.pt();
1338     jets_.genpeta[jets_.ngenp] = genp.eta();
1339     jets_.genpphi[jets_.ngenp] = genp.phi();
1340     jets_.genppdgId[jets_.ngenp] = genp.pdgId();
1341    
1342     jets_.ngenp++;
1343     }
1344     }
1345     }
1346    
1347    
1348     t->Fill();
1349    
1350    
1351    
1352     jets_.nj1 = 0;
1353     jets_.nj2 = 0;
1354     jets_.nj3 = 0;
1355     jets_.nj4 = 0;
1356     jets_.nPFcand = 0;
1357     jets_.ntrack = 0;
1358     jets_.ngenp = 0;
1359    
1360     }
1361    
1362    
1363     // copied from PhysicsTools/JetMCAlgos/plugins/PartonSelector.cc
1364     void HiPFJetAnalyzer::getPartons( const Event& iEvent, const EventSetup& iEs )
1365     {
1366    
1367     //edm::Handle <reco::CandidateView> genParticles;
1368     edm::Handle <reco::GenParticleCollection> genParticles;
1369     iEvent.getByLabel (genParticleTag_, genParticles );
1370    
1371     auto_ptr<GenParticleRefVector> thePartons ( new GenParticleRefVector);
1372    
1373    
1374     const GenParticle & parton1 = (*genParticles)[ 6 ];
1375     jets_.parton1_flavor = abs(parton1.pdgId());
1376     jets_.parton1_pt = parton1.pt();
1377     jets_.parton1_phi = parton1.phi();
1378     jets_.parton1_eta = parton1.eta();
1379     jets_.parton1_y = parton1.y();
1380    
1381     const GenParticle & parton2 = (*genParticles)[ 7 ];
1382     jets_.parton2_flavor = abs(parton2.pdgId());
1383     jets_.parton2_pt = parton2.pt();
1384     jets_.parton2_phi = parton2.phi();
1385     jets_.parton2_eta = parton2.eta();
1386     jets_.parton2_y = parton2.y();
1387    
1388    
1389    
1390    
1391     }
1392    
1393    
1394     //--------------------------------------------------------------------------------------------------
1395     void HiPFJetAnalyzer::fillL1Bits(const edm::Event &iEvent)
1396     {
1397     edm::Handle< L1GlobalTriggerReadoutRecord > L1GlobalTrigger;
1398    
1399     iEvent.getByLabel(L1gtReadout_, L1GlobalTrigger);
1400     const TechnicalTriggerWord& technicalTriggerWordBeforeMask = L1GlobalTrigger->technicalTriggerWord();
1401    
1402     for (int i=0; i<64;i++)
1403     {
1404     jets_.l1TBit[i] = technicalTriggerWordBeforeMask.at(i);
1405     }
1406     jets_.nL1TBit = 64;
1407    
1408     int ntrigs = L1GlobalTrigger->decisionWord().size();
1409     jets_.nL1ABit = ntrigs;
1410    
1411     for (int i=0; i != ntrigs; i++) {
1412     bool accept = L1GlobalTrigger->decisionWord()[i];
1413     //jets_.l1ABit[i] = (accept == true)? 1:0;
1414     if(accept== true){
1415     jets_.l1ABit[i] = 1;
1416     }
1417     else{
1418     jets_.l1ABit[i] = 0;
1419     }
1420    
1421     }
1422     }
1423    
1424     //--------------------------------------------------------------------------------------------------
1425     void HiPFJetAnalyzer::fillHLTBits(const edm::Event &iEvent)
1426     {
1427     // Fill HLT trigger bits.
1428     Handle<TriggerResults> triggerResultsHLT;
1429     getProduct(hltResName_, triggerResultsHLT, iEvent);
1430    
1431     const TriggerResults *hltResults = triggerResultsHLT.product();
1432     const TriggerNames & triggerNames = iEvent.triggerNames(*hltResults);
1433    
1434     jets_.nHLTBit = triggerNames.size();
1435    
1436     for(size_t i=0;i<hltTrgNames_.size();i++){
1437    
1438     for(size_t j=0;j<triggerNames.size();++j) {
1439    
1440     if(triggerNames.triggerName(j) == hltTrgNames_[i]){
1441    
1442     //cout <<"hltTrgNames_(i) "<<hltTrgNames_[i]<<endl;
1443     //cout <<"triggerName(j) "<<triggerNames.triggerName(j)<<endl;
1444     //cout<<" result "<<triggerResultsHLT->accept(j)<<endl;
1445     jets_.hltBit[i] = triggerResultsHLT->accept(j);
1446     }
1447    
1448     }
1449     }
1450     }
1451    
1452     //--------------------------------------------------------------------------------------------------
1453     template <typename TYPE>
1454     inline void HiPFJetAnalyzer::getProduct(const std::string name, edm::Handle<TYPE> &prod,
1455     const edm::Event &event) const
1456     {
1457     // Try to access data collection from EDM file. We check if we really get just one
1458     // product with the given name. If not we throw an exception.
1459    
1460     event.getByLabel(edm::InputTag(name),prod);
1461     if (!prod.isValid())
1462     throw edm::Exception(edm::errors::Configuration, "HiPFJetAnalyzer::GetProduct()\n")
1463     << "Collection with label '" << name << "' is not valid" << std::endl;
1464     }
1465    
1466     //--------------------------------------------------------------------------------------------------
1467     template <typename TYPE>
1468     inline bool HiPFJetAnalyzer::getProductSafe(const std::string name, edm::Handle<TYPE> &prod,
1469     const edm::Event &event) const
1470     {
1471     // Try to safely access data collection from EDM file. We check if we really get just one
1472     // product with the given name. If not, we return false.
1473    
1474     if (name.size()==0)
1475     return false;
1476    
1477     try {
1478     event.getByLabel(edm::InputTag(name),prod);
1479     if (!prod.isValid())
1480     return false;
1481     } catch (...) {
1482     return false;
1483     }
1484     return true;
1485     }
1486    
1487    
1488    
1489    
1490    
1491    
1492     DEFINE_FWK_MODULE(HiPFJetAnalyzer);