CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFJetTesterUnCorr.cc
Go to the documentation of this file.
6 
8 
11 
15 
16 //I don't know which of these I actually need yet
23 
25 
27 
28 #include "PFJetTesterUnCorr.h"
29 
30 #include <cmath>
31 
32 using namespace edm;
33 using namespace reco;
34 using namespace std;
35 
36 namespace {
37  bool is_B (const reco::Jet& fJet) {return fabs (fJet.eta()) < 1.3;}
38  bool is_E (const reco::Jet& fJet) {return fabs (fJet.eta()) >= 1.3 && fabs (fJet.eta()) < 3.;}
39  bool is_F (const reco::Jet& fJet) {return fabs (fJet.eta()) >= 3.;}
40 }
41 
43  : mInputCollection (iConfig.getParameter<edm::InputTag>( "src" )),
44  mInputGenCollection (iConfig.getParameter<edm::InputTag>( "srcGen" )),
45  rho_tag_ (iConfig.getParameter<edm::InputTag>( "srcRho" )),
46  mOutputFile (iConfig.getUntrackedParameter<std::string>("outputFile", "")),
47  mMatchGenPtThreshold (iConfig.getParameter<double>("genPtThreshold")),
48  mGenEnergyFractionThreshold (iConfig.getParameter<double>("genEnergyFractionThreshold")),
49  mReverseEnergyFractionThreshold (iConfig.getParameter<double>("reverseEnergyFractionThreshold")),
50  mRThreshold (iConfig.getParameter<double>("RThreshold")),
51  mTurnOnEverything (iConfig.getUntrackedParameter<std::string>("TurnOnEverything",""))
52 {
55  = mP = mP_80 = mPt = mPt_80
59  // = mMaxEInEmTowers = mMaxEInHadTowers
60  // = mHadEnergyInHO = mHadEnergyInHB = mHadEnergyInHE
62  // = mHadEnergyInHO_80 = mHadEnergyInHB_80 = mHadEnergyInHE_80
63  // = mHadEnergyInHO_3000 = mHadEnergyInHB_3000 = mHadEnergyInHE_3000
64  // = mEmEnergyInEB = mEmEnergyInEE
66  // = mEmEnergyInEB_80 = mEmEnergyInEE_80
67  // = mEmEnergyInEB_3000 = mEmEnergyInEE_3000
68  // = mEnergyFractionHadronic = mEnergyFractionEm
69  // = mHFLong = mHFTotal = mHFLong_80 = mHFLong_3000 = mHFShort = mHFShort_80 = mHFShort_3000
70  // = mN90
77  //= mCaloMEx = mCaloMEx_3000 = mCaloMEy = mCaloMEy_3000 = mCaloMETSig = mCaloMETSig_3000
78  //= mCaloMET = mCaloMET_3000 = mCaloMETPhi = mCaloSumET = mCaloSumET_3000
81  //= mAllGenJetsPt = mMatchedGenJetsPt = mAllGenJetsEta = mMatchedGenJetsEta
82  //= mGenJetMatchEnergyFraction = mReverseMatchEnergyFraction = mRMatch
83  = mDeltaEta = mDeltaPhi //= mEScale = mlinEScale = mDeltaE
84  //= mHadEnergyProfile = mEmEnergyProfile = mJetEnergyProfile
85  // = mHadJetEnergyProfile = mEMJetEnergyProfile
87  //= mpTScaleB_s = mpTScaleE_s = mpTScaleF_s
90  //= mpTScale_30_200_s = mpTScale_200_600_s = mpTScale_600_1500_s = mpTScale_1500_3500_s
92 
97  /*
98  = mpTScale1D_30_200 = mpTScale1D_200_600 = mpTScale1D_600_1500 = mpTScale1D_1500_3500
99  = mHBEne = mHBTime = mHEEne = mHETime = mHFEne = mHFTime = mHOEne = mHOTime
100  = mEBEne = mEBTime = mEEEne = mEETime
101  */
102 
104  = 0;
105 
106  DQMStore* dbe = &*edm::Service<DQMStore>();
107  if (dbe) {
108  dbe->setCurrentFolder("JetMET/RecoJetsV/PFJetTask_" + mInputCollection.label());
109  //
110  numberofevents = dbe->book1D("numberofevents","numberofevents", 3, 0 , 2);
111  //
112  mEta = dbe->book1D("Eta", "Eta", 120, -6, 6);
113  mEtaFineBin = dbe->book1D("EtaFineBin_Pt10", "EtaFineBin_Pt10", 600, -6, 6);
114  /*
115  mEtaFineBin1p = dbe->book1D("EtaFineBin1p_Pt10", "EtaFineBin1p_Pt10", 100, 0, 1.3);
116  mEtaFineBin2p = dbe->book1D("EtaFineBin2p_Pt10", "EtaFineBin2p_Pt10", 100, 1.3, 3);
117  mEtaFineBin3p = dbe->book1D("EtaFineBin3p_Pt10", "EtaFineBin3p_Pt10", 100, 3, 5);
118  mEtaFineBin1m = dbe->book1D("EtaFineBin1m_Pt10", "EtaFineBin1m_Pt10", 100, -1.3, 0);
119  mEtaFineBin2m = dbe->book1D("EtaFineBin2m_Pt10", "EtaFineBin2m_Pt10", 100, -3, -1.3);
120  mEtaFineBin3m = dbe->book1D("EtaFineBin3m_Pt10", "EtaFineBin3m_Pt10", 100, -5, -3);
121  */
122  //
123  mPhi = dbe->book1D("Phi", "Phi", 70, -3.5, 3.5);
124  mPhiFineBin = dbe->book1D("PhiFineBin_Pt10", "PhiFineBin_Pt10", 350, -3.5, 3.5);
125  //
126  mE = dbe->book1D("E", "E", 100, 0, 500);
127  mE_80 = dbe->book1D("E_80", "E_80", 100, 0, 5000);
128  //
129  mP = dbe->book1D("P", "P", 100, 0, 500);
130  mP_80 = dbe->book1D("P_80", "P_80", 100, 0, 5000);
131  //
132  mPt = dbe->book1D("Pt", "Pt", 100, 0, 150);
133  mPt_80 = dbe->book1D("Pt_80", "Pt_80", 100, 0, 4000);
134  //
135  mMass = dbe->book1D("Mass", "Mass", 100, 0, 200);
136  mMass_80 = dbe->book1D("Mass_80", "Mass_80", 100, 0, 500);
137  //
138  mConstituents = dbe->book1D("Constituents", "# of Constituents", 100, 0, 100);
139  mConstituents_80 = dbe->book1D("Constituents_80", "# of Constituents_80", 40, 0, 40);
140  //
141  mEtaFirst = dbe->book1D("EtaFirst", "EtaFirst", 120, -6, 6);
142  mPhiFirst = dbe->book1D("PhiFirst", "PhiFirst", 70, -3.5, 3.5);
143  mPtFirst = dbe->book1D("PtFirst", "PtFirst", 100, 0, 50);
144  mPtFirst_80 = dbe->book1D("PtFirst_80", "PtFirst_80", 100, 0, 140);
145  mPtFirst_3000 = dbe->book1D("PtFirst_3000", "PtFirst_3000", 100, 0, 4000);
146  //
147  mMjj = dbe->book1D("Mjj", "Mjj", 100, 0, 2000);
148  mMjj_3000 = dbe->book1D("Mjj_3000", "Mjj_3000", 100, 0, 10000);
149  mDelEta = dbe->book1D("DelEta", "DelEta", 100, -.5, .5);
150  mDelPhi = dbe->book1D("DelPhi", "DelPhi", 100, -.5, .5);
151  mDelPt = dbe->book1D("DelPt", "DelPt", 100, -1, 1);
152  //
153  // mMaxEInEmTowers = dbe->book1D("MaxEInEmTowers", "MaxEInEmTowers", 100, 0, 100);
154  // mMaxEInHadTowers = dbe->book1D("MaxEInHadTowers", "MaxEInHadTowers", 100, 0, 100);
155  // mHadEnergyInHO = dbe->book1D("HadEnergyInHO", "HadEnergyInHO", 100, 0, 10);
156  // mHadEnergyInHB = dbe->book1D("HadEnergyInHB", "HadEnergyInHB", 100, 0, 50);
157  mHadEnergyInHF = dbe->book1D("HadEnergyInHF", "HadEnergyInHF", 100, 0, 2500);
158  mHadEnergyInHF_80 = dbe->book1D("HadEnergyInHF_80", "HadEnergyInHF_80", 100, 0, 3000);
159  mHadEnergyInHF_3000 = dbe->book1D("HadEnergyInHF_3000", "HadEnergyInHF_3000", 100, 0, 1800);
160  // mHadEnergyInHE = dbe->book1D("HadEnergyInHE", "HadEnergyInHE", 100, 0, 100);
161  //
162  // mHadEnergyInHO_80 = dbe->book1D("HadEnergyInHO_80", "HadEnergyInHO_80", 100, 0, 50);
163  // mHadEnergyInHB_80 = dbe->book1D("HadEnergyInHB_80", "HadEnergyInHB_80", 100, 0, 200);
164  // mHadEnergyInHE_80 = dbe->book1D("HadEnergyInHE_80", "HadEnergyInHE_80", 100, 0, 1000);
165  // mHadEnergyInHO_3000 = dbe->book1D("HadEnergyInHO_3000", "HadEnergyInHO_3000", 100, 0, 500);
166  // mHadEnergyInHB_3000 = dbe->book1D("HadEnergyInHB_3000", "HadEnergyInHB_3000", 100, 0, 3000);
167  // mHadEnergyInHE_3000 = dbe->book1D("HadEnergyInHE_3000", "HadEnergyInHE_3000", 100, 0, 2000);
168  //
169  // mEmEnergyInEB = dbe->book1D("EmEnergyInEB", "EmEnergyInEB", 100, 0, 50);
170  // mEmEnergyInEE = dbe->book1D("EmEnergyInEE", "EmEnergyInEE", 100, 0, 50);
171  mEmEnergyInHF = dbe->book1D("EmEnergyInHF", "EmEnergyInHF", 100, -20, 450);
172  mEmEnergyInHF_80 = dbe->book1D("EmEnergyInHF_80", "EmEnergyInHF_80", 100, -20, 440);
173  mEmEnergyInHF_3000 = dbe->book1D("EmEnergyInHF_3000", "EmEnergyInHF_3000", 100, -20, 190);
174  // mEmEnergyInEB_80 = dbe->book1D("EmEnergyInEB_80", "EmEnergyInEB_80", 100, 0, 200);
175  // mEmEnergyInEE_80 = dbe->book1D("EmEnergyInEE_80", "EmEnergyInEE_80", 100, 0, 1000);
176  // mEmEnergyInEB_3000= dbe->book1D("EmEnergyInEB_3000", "EmEnergyInEB_3000", 100, 0, 3000);
177  // mEmEnergyInEE_3000= dbe->book1D("EmEnergyInEE_3000", "EmEnergyInEE_3000", 100, 0, 2000);
178  // mEnergyFractionHadronic = dbe->book1D("EnergyFractionHadronic", "EnergyFractionHadronic", 120, -0.1, 1.1);
179  // mEnergyFractionEm = dbe->book1D("EnergyFractionEm", "EnergyFractionEm", 120, -0.1, 1.1);
180  //
181  // mHFTotal = dbe->book1D("HFTotal", "HFTotal", 100, 0, 500);
182  // mHFTotal_80 = dbe->book1D("HFTotal_80", "HFTotal_80", 100, 0, 3000);
183  // mHFTotal_3000 = dbe->book1D("HFTotal_3000", "HFTotal_3000", 100, 0, 6000);
184  // mHFLong = dbe->book1D("HFLong", "HFLong", 100, 0, 500);
185  // mHFLong_80 = dbe->book1D("HFLong_80", "HFLong_80", 100, 0, 200);
186  // mHFLong_3000 = dbe->book1D("HFLong_3000", "HFLong_3000", 100, 0, 1500);
187  // mHFShort = dbe->book1D("HFShort", "HFShort", 100, 0, 500);
188  // mHFShort_80 = dbe->book1D("HFShort_80", "HFShort_80", 100, 0, 200);
189  // mHFShort_3000 = dbe->book1D("HFShort_3000", "HFShort_3000", 100, 0, 1500);
190  //
191  // mN90 = dbe->book1D("N90", "N90", 50, 0, 50);
192  //
193 
194  mChargedEmEnergy_80 = dbe->book1D("ChargedEmEnergy_80","ChargedEmEnergy_80",100,0.,500.);
195  mChargedHadronEnergy_80 = dbe->book1D("ChargedHadronEnergy_80","ChargedHadronEnergy_80",100,0.,500.);
196  mNeutralEmEnergy_80 = dbe->book1D("NeutralEmEnergy_80","NeutralEmEnergy_80",100,0.,500.);
197  mNeutralHadronEnergy_80 = dbe->book1D("NeutralHadronEnergy_80","NeutralHadronEnergy_80",100,0.,500.);
198 
199  mChargedEmEnergy_3000 = dbe->book1D("ChargedEmEnergy_3000","ChargedEmEnergy_3000",100,0.,2000.);
200  mChargedHadronEnergy_3000 = dbe->book1D("ChargedHadronEnergy_3000","ChargedHadronEnergy_3000",100,0.,2000.);
201  mNeutralEmEnergy_3000 = dbe->book1D("NeutralEmEnergy_3000","NeutralEmEnergy_3000",100,0.,2000.);
202  mNeutralHadronEnergy_3000 = dbe->book1D("NeutralHadronEnergy_3000","NeutralHadronEnergy_3000",100,0.,2000.);
203 
204  mChargedEmEnergyFraction_B = dbe->book1D("ChargedEmEnergyFraction_B","ChargedEmEnergyFraction_B",120,-0.1,1.1);
205  mChargedEmEnergyFraction_E = dbe->book1D("ChargedEmEnergyFraction_E","ChargedEmEnergyFraction_E",120,-0.1,1.1);
206  mChargedEmEnergyFraction_F = dbe->book1D("ChargedEmEnergyFraction_F","ChargedEmEnergyFraction_F",120,-0.1,1.1);
207  mChargedHadronEnergyFraction_B = dbe->book1D("ChargedHadronEnergyFraction_B","ChargedHadronEnergyFraction_B",120,-0.1,1.1);
208  mChargedHadronEnergyFraction_E = dbe->book1D("ChargedHadronEnergyFraction_E","ChargedHadronEnergyFraction_E",120,-0.1,1.1);
209  mChargedHadronEnergyFraction_F = dbe->book1D("ChargedHadronEnergyFraction_F","ChargedHadronEnergyFraction_F",120,-0.1,1.1);
210  mNeutralEmEnergyFraction_B = dbe->book1D("NeutralEmEnergyFraction_B","NeutralEmEnergyFraction_B",120,-0.1,1.1);
211  mNeutralEmEnergyFraction_E = dbe->book1D("NeutralEmEnergyFraction_E","NeutralEmEnergyFraction_E",120,-0.1,1.1);
212  mNeutralEmEnergyFraction_F = dbe->book1D("NeutralEmEnergyFraction_F","NeutralEmEnergyFraction_F",120,-0.1,1.1);
213  mNeutralHadronEnergyFraction_B = dbe->book1D("NeutralHadronEnergyFraction_B","NeutralHadronEnergyFraction_B",120,-0.1,1.1);
214  mNeutralHadronEnergyFraction_E = dbe->book1D("NeutralHadronEnergyFraction_E","NeutralHadronEnergyFraction_E",120,-0.1,1.1);
215  mNeutralHadronEnergyFraction_F = dbe->book1D("NeutralHadronEnergyFraction_F","NeutralHadronEnergyFraction_F",120,-0.1,1.1);
216 
217  mGenEta = dbe->book1D("GenEta", "GenEta", 120, -6, 6);
218  mGenPhi = dbe->book1D("GenPhi", "GenPhi", 70, -3.5, 3.5);
219  mGenPt = dbe->book1D("GenPt", "GenPt", 100, 0, 150);
220  mGenPt_80 = dbe->book1D("GenPt_80", "GenPt_80", 100, 0, 1500);
221  //
222  mGenEtaFirst = dbe->book1D("GenEtaFirst", "GenEtaFirst", 120, -6, 6);
223  mGenPhiFirst = dbe->book1D("GenPhiFirst", "GenPhiFirst", 70, -3.5, 3.5);
224  //
225  /*
226  mCaloMEx = dbe->book1D("CaloMEx","CaloMEx",200,-150,150);
227  mCaloMEx_3000 = dbe->book1D("CaloMEx_3000","CaloMEx_3000",100,-500,500);
228  mCaloMEy = dbe->book1D("CaloMEy","CaloMEy",200,-150,150);
229  mCaloMEy_3000 = dbe->book1D("CaloMEy_3000","CaloMEy_3000",100,-500,500);
230  mCaloMETSig = dbe->book1D("CaloMETSig","CaloMETSig",100,0,15);
231  mCaloMETSig_3000 = dbe->book1D("CaloMETSig_3000","CaloMETSig_3000",100,0,50);
232  mCaloMET = dbe->book1D("CaloMET","CaloMET",100,0,150);
233  mCaloMET_3000 = dbe->book1D("CaloMET_3000","CaloMET_3000",100,0,1000);
234  mCaloMETPhi = dbe->book1D("CaloMETPhi","CaloMETPhi",70, -3.5, 3.5);
235  mCaloSumET = dbe->book1D("CaloSumET","CaloSumET",100,0,500);
236  mCaloSumET_3000 = dbe->book1D("CaloSumET_3000","CaloSumET_3000",100,3000,8000);
237  */
238  //
239  mHadTiming = dbe->book1D("HadTiming", "HadTiming", 75, -50, 100);
240  mEmTiming = dbe->book1D("EMTiming", "EMTiming", 75, -50, 100);
241  //
242  mNJetsEtaC = dbe->book1D("NJetsEtaC_Pt10", "NJetsEtaC_Pt10", 15, 0, 15);
243  mNJetsEtaF = dbe->book1D("NJetsEtaF_Pt10", "NJetsEtaF_Pt10", 15, 0, 15);
244  //
245  mNJets1 = dbe->bookProfile("NJets1", "NJets1", 100, 0, 200, 100, 0, 50, "s");
246  mNJets2 = dbe->bookProfile("NJets2", "NJets2", 100, 0, 4000, 100, 0, 50, "s");
247  //
248  /*
249  mHBEne = dbe->book1D( "HBEne", "HBEne", 1000, -20, 100 );
250  mHBTime = dbe->book1D( "HBTime", "HBTime", 200, -200, 200 );
251  mHEEne = dbe->book1D( "HEEne", "HEEne", 1000, -20, 100 );
252  mHETime = dbe->book1D( "HETime", "HETime", 200, -200, 200 );
253  mHOEne = dbe->book1D( "HOEne", "HOEne", 1000, -20, 100 );
254  mHOTime = dbe->book1D( "HOTime", "HOTime", 200, -200, 200 );
255  mHFEne = dbe->book1D( "HFEne", "HFEne", 1000, -20, 100 );
256  mHFTime = dbe->book1D( "HFTime", "HFTime", 200, -200, 200 );
257  mEBEne = dbe->book1D( "EBEne", "EBEne", 1000, -20, 100 );
258  mEBTime = dbe->book1D( "EBTime", "EBTime", 200, -200, 200 );
259  mEEEne = dbe->book1D( "EEEne", "EEEne", 1000, -20, 100 );
260  mEETime = dbe->book1D( "EETime", "EETime", 200, -200, 200 );
261  */
262  //
263  mPthat_80 = dbe->book1D("Pthat_80", "Pthat_80", 100, 0.0, 1000.0);
264  mPthat_3000 = dbe->book1D("Pthat_3000", "Pthat_3000", 100, 1000.0, 4000.0);
265 
266  mjetArea = dbe->book1D("jetArea","jetArea",26,-0.5,12.5);
267  mRho = dbe->book1D("Rho","Rho",20,0,20);
268 
269  //
270  double log10PtMin = 0.5; //=3.1622766
271  double log10PtMax = 3.75; //=5623.41325
272  int log10PtBins = 26;
273  //double etaMin = -6;
274  //double etaMax = 6.;
275  //int etaBins = 50;
276  double etaRange[91] = {-6.0,-5.8,-5.6,-5.4,-5.2,-5.0,-4.8,-4.6,-4.4,-4.2,-4.0,-3.8,-3.6,-3.4,-3.2,-3.0,-2.9,-2.8,-2.7,-2.6,-2.5,-2.4,-2.3,-2.2,-2.1,-2.0,-1.9,-1.8,-1.7,-1.6,-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,3.0,3.2,3.4,3.6,3.8,4.0,4.2,4.4,4.6,4.8,5.0,5.2,5.4,5.6,5.8,6.0};
277  //double linPtMin = 5;
278  //double linPtMax = 155;
279  //int linPtBins = 15;
280 
281  //int log10PtFineBins = 50;
282  /*
283  mAllGenJetsPt = dbe->book1D("GenJetLOGpT", "GenJet LOG(pT_gen)",
284  log10PtBins, log10PtMin, log10PtMax);
285  mMatchedGenJetsPt = dbe->book1D("MatchedGenJetLOGpT", "MatchedGenJet LOG(pT_gen)",
286  log10PtBins, log10PtMin, log10PtMax);
287  mAllGenJetsEta = dbe->book2D("GenJetEta", "GenJet Eta vs LOG(pT_gen)",
288  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax);
289  mMatchedGenJetsEta = dbe->book2D("MatchedGenJetEta", "MatchedGenJet Eta vs LOG(pT_gen)",
290  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax);
291  */
292  //
293  if (mTurnOnEverything.compare("yes")==0) {
294  /*
295  mHadEnergyProfile = dbe->bookProfile2D("HadEnergyProfile", "HadEnergyProfile", 82, -41, 41, 73, 0, 73, 100, 0, 10000, "s");
296  mEmEnergyProfile = dbe->bookProfile2D("EmEnergyProfile", "EmEnergyProfile", 82, -41, 41, 73, 0, 73, 100, 0, 10000, "s");
297  */
298  }
299  //mJetEnergyProfile = dbe->bookProfile2D("JetEnergyProfile", "JetEnergyProfile", 50, -5, 5, 36, -3.1415987, 3.1415987, 100, 0, 10000, "s");
300  // mHadJetEnergyProfile = dbe->bookProfile2D("HadJetEnergyProfile", "HadJetEnergyProfile", 50, -5, 5, 36, -3.1415987, 3.1415987, 100, 0, 10000, "s");
301  // mEMJetEnergyProfile = dbe->bookProfile2D("EMJetEnergyProfile", "EMJetEnergyProfile", 50, -5, 5, 36, -3.1415987, 3.1415987, 100, 0, 10000, "s");
302  //
303  if (mTurnOnEverything.compare("yes")==0) {
304  /*
305  mGenJetMatchEnergyFraction = dbe->book3D("GenJetMatchEnergyFraction", "GenJetMatchEnergyFraction vs LOG(pT_gen) vs eta",
306  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 101, 0, 1.01);
307  mReverseMatchEnergyFraction = dbe->book3D("ReverseMatchEnergyFraction", "ReverseMatchEnergyFraction vs LOG(pT_gen) vs eta",
308  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 101, 0, 1.01);
309  mRMatch = dbe->book3D("RMatch", "delta(R)(Gen-Calo) vs LOG(pT_gen) vs eta",
310  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 60, 0, 3);
311  */
312 /*
313  mDeltaEta = dbe->book3D("DeltaEta", "DeltaEta vs LOG(pT_gen) vs eta",
314  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 100, -1, 1);
315  mDeltaPhi = dbe->book3D("DeltaPhi", "DeltaPhi vs LOG(pT_gen) vs eta",
316  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 100, -1, 1);
317 */
318  /*
319  mEScale = dbe->book3D("EScale", "EnergyScale vs LOG(pT_gen) vs eta",
320  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 100, 0, 2);
321  mlinEScale = dbe->book3D("linEScale", "EnergyScale vs LOG(pT_gen) vs eta",
322  linPtBins, linPtMin, linPtMax, etaBins, etaMin, etaMax, 100, 0, 2);
323  mDeltaE = dbe->book3D("DeltaE", "DeltaE vs LOG(pT_gen) vs eta",
324  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 2000, -200, 200);
325  */
326  //
327  /* mEScale_pt10 = dbe->book3D("EScale_pt10", "EnergyScale vs LOG(pT_gen) vs eta",
328  log10PtBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 100, 0, 2);
329  mEScaleFineBin = dbe->book3D("EScaleFineBins", "EnergyScale vs LOG(pT_gen) vs eta",
330  log10PtFineBins, log10PtMin, log10PtMax, etaBins, etaMin, etaMax, 100, 0, 2);
331 */
332  }
333  /*
334  mpTScaleB_s = dbe->bookProfile("pTScaleB_s", "pTScale_s_0<|eta|<1.3",
335  log10PtBins, log10PtMin, log10PtMax, 0, 2, "s");
336  mpTScaleE_s = dbe->bookProfile("pTScaleE_s", "pTScale_s_1.3<|eta|<3.0",
337  log10PtBins, log10PtMin, log10PtMax, 0, 2, "s");
338  mpTScaleF_s = dbe->bookProfile("pTScaleF_s", "pTScale_s_3.0<|eta|<5.0",
339  log10PtBins, log10PtMin, log10PtMax, 0, 2, "s");
340  */
341  mpTScaleB_d = dbe->bookProfile("pTScaleB_d", "pTScale_d_0<|eta|<1.5",
342  log10PtBins, log10PtMin, log10PtMax, 0, 2, " ");
343  mpTScaleE_d = dbe->bookProfile("pTScaleE_d", "pTScale_d_1.5<|eta|<3.0",
344  log10PtBins, log10PtMin, log10PtMax, 0, 2, " ");
345  mpTScaleF_d = dbe->bookProfile("pTScaleF_d", "pTScale_d_3.0<|eta|<6.0",
346  log10PtBins, log10PtMin, log10PtMax, 0, 2, " ");
347  mpTScalePhiB_d = dbe->bookProfile("pTScalePhiB_d", "pTScalePhi_d_0<|eta|<1.5",
348  70, -3.5, 3.5, 0, 2, " ");
349  mpTScalePhiE_d = dbe->bookProfile("pTScalePhiE_d", "pTScalePhi_d_1.5<|eta|<3.0",
350  70, -3.5, 3.5, 0, 2, " ");
351  mpTScalePhiF_d = dbe->bookProfile("pTScalePhiF_d", "pTScalePhi_d_3.0<|eta|<6.0",
352  70, -3.5, 3.5, 0, 2, " ");
353  /*
354  mpTScale_30_200_s = dbe->bookProfile("pTScale_30_200_s", "pTScale_s_30<pT<200",
355  etaBins, etaMin, etaMax, 0., 2., "s");
356  mpTScale_200_600_s = dbe->bookProfile("pTScale_200_600_s", "pTScale_s_200<pT<600",
357  etaBins, etaMin, etaMax, 0., 2., "s");
358  mpTScale_600_1500_s = dbe->bookProfile("pTScale_600_1500_s", "pTScale_s_600<pT<1500",
359  etaBins, etaMin, etaMax, 0., 2., "s");
360  mpTScale_1500_3500_s = dbe->bookProfile("pTScale_1500_3500_s", "pTScale_s_1500<pt<3500",
361  etaBins, etaMin, etaMax, 0., 2., "s");
362  */
363  mpTScale_30_200_d = dbe->bookProfile("pTScale_30_200_d", "pTScale_d_30<pT<200",
364  90,etaRange, 0., 2., " ");
365  mpTScale_200_600_d = dbe->bookProfile("pTScale_200_600_d", "pTScale_d_200<pT<600",
366  90,etaRange, 0., 2., " ");
367  mpTScale_600_1500_d = dbe->bookProfile("pTScale_600_1500_d", "pTScale_d_600<pT<1500",
368  90,etaRange, 0., 2., " ");
369  mpTScale_1500_3500_d = dbe->bookProfile("pTScale_1500_3500_d", "pTScale_d_1500<pt<3500",
370  90,etaRange, 0., 2., " ");
371 
372  mpTScale1DB_30_200 = dbe->book1D("pTScale1DB_30_200", "pTScale_distribution_for_0<|eta|<1.5_30_200",
373  100, 0, 2);
374  mpTScale1DE_30_200 = dbe->book1D("pTScale1DE_30_200", "pTScale_distribution_for_1.5<|eta|<3.0_30_200",
375  50, 0, 2);
376  mpTScale1DF_30_200 = dbe->book1D("pTScale1DF_30_200", "pTScale_distribution_for_3.0<|eta|<6.0_30_200",
377  50, 0, 2);
378 
379  mpTScale1DB_200_600 = dbe->book1D("pTScale1DB_200_600", "pTScale_distribution_for_0<|eta|<1.5_200_600",
380  100, 0, 2);
381  mpTScale1DE_200_600 = dbe->book1D("pTScale1DE_200_600", "pTScale_distribution_for_1.5<|eta|<3.0_200_600",
382  50, 0, 2);
383  mpTScale1DF_200_600 = dbe->book1D("pTScale1DF_200_600", "pTScale_distribution_for_3.0<|eta|<6.0_200_600",
384  50, 0, 2);
385 
386  mpTScale1DB_600_1500 = dbe->book1D("pTScale1DB_600_1500", "pTScale_distribution_for_0<|eta|<1.5_600_1500",
387  100, 0, 2);
388  mpTScale1DE_600_1500 = dbe->book1D("pTScale1DE_600_1500", "pTScale_distribution_for_1.5<|eta|<3.0_600_1500",
389  50, 0, 2);
390  mpTScale1DF_600_1500 = dbe->book1D("pTScale1DF_600_1500", "pTScale_distribution_for_3.0<|eta|<6.0_600_1500",
391  50, 0, 2);
392 
393  mpTScale1DB_1500_3500 = dbe->book1D("pTScale1DB_1500_3500", "pTScale_distribution_for_0<|eta|<1.5_1500_3500",
394  100, 0, 2);
395  mpTScale1DE_1500_3500 = dbe->book1D("pTScale1DE_1500_3500", "pTScale_distribution_for_1.5<|eta|<3.0_1500_3500",
396  50, 0, 2);
397  mpTScale1DF_1500_3500 = dbe->book1D("pTScale1DF_1500_3500", "pTScale_distribution_for_3.0<|eta|<6.0_1500_3500",
398  50, 0, 2);
399  /*
400  mpTScale1D_30_120 = dbe->book1D("pTScale1D_30_120", "pTScale_distribution_for_30<pT<120",
401  100, 0, 2);
402  mpTScale1D_200_600 = dbe->book1D("pTScale1D_200_600", "pTScale_distribution_for_200<pT<600",
403  100, 0, 2);
404  mpTScale1D_600_1500 = dbe->book1D("pTScale1D_600_1500", "pTScale_distribution_for_600<pT<1500",
405  100, 0, 2);
406  mpTScale1D_1500_3500 = dbe->book1D("pTScale1D_1500_3500", "pTScale_distribution_for_1500<pt<3500",
407  100, 0, 2);
408  */
409 
410  }
411 
412  if (mOutputFile.empty ()) {
413  LogInfo("OutputInfo") << " PFJet histograms will NOT be saved";
414  }
415  else {
416  LogInfo("OutputInfo") << " PFJethistograms will be saved to file:" << mOutputFile;
417  }
418 }
419 
421 {
422 }
423 
425 }
426 
429 }
430 
431 
432 void PFJetTesterUnCorr::analyze(const edm::Event& mEvent, const edm::EventSetup& mSetup)
433 {
434  double countsfornumberofevents = 1;
435  numberofevents->Fill(countsfornumberofevents);
436  // *********************************
437  // *** Get pThat
438  // *********************************
439 if (!mEvent.isRealData()){
441  mEvent.getByLabel("generator", evt);
442  if (evt.isValid()) {
443  HepMC::GenEvent * myGenEvent = new HepMC::GenEvent(*(evt->GetEvent()));
444 
445  double pthat = myGenEvent->event_scale();
446 
447  mPthat_80->Fill(pthat);
448  mPthat_3000->Fill(pthat);
449 
450  delete myGenEvent;
451  }
452 }
453  // ***********************************
454  // *** Get CaloMET
455  // ***********************************
456 /*
457  const CaloMET *calomet;
458  edm::Handle<CaloMETCollection> calo;
459  mEvent.getByLabel("met", calo);
460  if (!calo.isValid()) {
461  edm::LogInfo("OutputInfo") << " failed to retrieve data required by MET Task";
462  edm::LogInfo("OutputInfo") << " MET Task cannot continue...!";
463  } else {
464  const CaloMETCollection *calometcol = calo.product();
465  calomet = &(calometcol->front());
466 
467  double caloSumET = calomet->sumEt();
468  double caloMETSig = calomet->mEtSig();
469  double caloMET = calomet->pt();
470  double caloMEx = calomet->px();
471  double caloMEy = calomet->py();
472  double caloMETPhi = calomet->phi();
473 
474  mCaloMEx->Fill(caloMEx);
475  mCaloMEx_3000->Fill(caloMEx);
476  mCaloMEy->Fill(caloMEy);
477  mCaloMEy_3000->Fill(caloMEy);
478  mCaloMET->Fill(caloMET);
479  mCaloMET_3000->Fill(caloMET);
480  mCaloMETPhi->Fill(caloMETPhi);
481  mCaloSumET->Fill(caloSumET);
482  mCaloSumET_3000->Fill(caloSumET);
483  mCaloMETSig->Fill(caloMETSig);
484  mCaloMETSig_3000->Fill(caloMETSig);
485 
486  }
487 */
488  // ***********************************
489  // *** Get the CaloTower collection
490  // ***********************************
492  mEvent.getByLabel( "towerMaker", caloTowers );
493  if (caloTowers.isValid()) {
494  for( CaloTowerCollection::const_iterator cal = caloTowers->begin(); cal != caloTowers->end(); ++ cal ){
495 
496  //To compensate for the index
497  if (mTurnOnEverything.compare("yes")==0) {
498  /*
499  if (cal->ieta() >> 0 ){mHadEnergyProfile->Fill (cal->ieta()-1, cal->iphi(), cal->hadEnergy());
500  mEmEnergyProfile->Fill (cal->ieta()-1, cal->iphi(), cal->emEnergy());}
501  mHadEnergyProfile->Fill (cal->ieta(), cal->iphi(), cal->hadEnergy());
502  mEmEnergyProfile->Fill (cal->ieta(), cal->iphi(), cal->emEnergy());
503  */
504  }
505 
506  mHadTiming->Fill (cal->hcalTime());
507  mEmTiming->Fill (cal->ecalTime());
508  }
509  }
510 
511  // ***********************************
512  // *** Get the RecHits collection
513  // ***********************************
514  try {
515  std::vector<edm::Handle<HBHERecHitCollection> > colls;
516  mEvent.getManyByType(colls);
517  std::vector<edm::Handle<HBHERecHitCollection> >::iterator i;
518  for (i=colls.begin(); i!=colls.end(); i++) {
519  for (HBHERecHitCollection::const_iterator j=(*i)->begin(); j!=(*i)->end(); j++) {
520  // std::cout << *j << std::endl;
521  /*
522  if (j->id().subdet() == HcalBarrel) {
523  mHBEne->Fill(j->energy());
524  mHBTime->Fill(j->time());
525  }
526  if (j->id().subdet() == HcalEndcap) {
527  mHEEne->Fill(j->energy());
528  mHETime->Fill(j->time());
529  }
530  */
531  }
532  }
533  } catch (...) {
534  edm::LogInfo("OutputInfo") << " No HB/HE RecHits.";
535  }
536 
537  try {
538  std::vector<edm::Handle<HFRecHitCollection> > colls;
539  mEvent.getManyByType(colls);
540  std::vector<edm::Handle<HFRecHitCollection> >::iterator i;
541  for (i=colls.begin(); i!=colls.end(); i++) {
542  for (HFRecHitCollection::const_iterator j=(*i)->begin(); j!=(*i)->end(); j++) {
543  // std::cout << *j << std::endl;
544  /*
545  if (j->id().subdet() == HcalForward) {
546  mHFEne->Fill(j->energy());
547  mHFTime->Fill(j->time());
548  }
549  */
550  }
551  }
552  } catch (...) {
553  edm::LogInfo("OutputInfo") << " No HF RecHits.";
554  }
555 
556  try {
557  std::vector<edm::Handle<HORecHitCollection> > colls;
558  mEvent.getManyByType(colls);
559  std::vector<edm::Handle<HORecHitCollection> >::iterator i;
560  for (i=colls.begin(); i!=colls.end(); i++) {
561  for (HORecHitCollection::const_iterator j=(*i)->begin(); j!=(*i)->end(); j++) {
562  /*
563  if (j->id().subdet() == HcalOuter) {
564  mHOEne->Fill(j->energy());
565  mHOTime->Fill(j->time());
566  }
567  */
568  }
569  }
570  } catch (...) {
571  edm::LogInfo("OutputInfo") << " No HO RecHits.";
572  }
573  try {
574  std::vector<edm::Handle<EBRecHitCollection> > colls;
575  mEvent.getManyByType(colls);
576  std::vector<edm::Handle<EBRecHitCollection> >::iterator i;
577  for (i=colls.begin(); i!=colls.end(); i++) {
578  for (EBRecHitCollection::const_iterator j=(*i)->begin(); j!=(*i)->end(); j++) {
579  // if (j->id() == EcalBarrel) {
580  //mEBEne->Fill(j->energy());
581  //mEBTime->Fill(j->time());
582  // }
583  // std::cout << *j << std::endl;
584  // std::cout << j->id() << std::endl;
585  }
586  }
587  } catch (...) {
588  edm::LogInfo("OutputInfo") << " No EB RecHits.";
589  }
590 
591  try {
592  std::vector<edm::Handle<EERecHitCollection> > colls;
593  mEvent.getManyByType(colls);
594  std::vector<edm::Handle<EERecHitCollection> >::iterator i;
595  for (i=colls.begin(); i!=colls.end(); i++) {
596  for (EERecHitCollection::const_iterator j=(*i)->begin(); j!=(*i)->end(); j++) {
597  // if (j->id().subdet() == EcalEndcap) {
598  //mEEEne->Fill(j->energy());
599  //mEETime->Fill(j->time());
600  // }
601  // std::cout << *j << std::endl;
602  }
603  }
604  } catch (...) {
605  edm::LogInfo("OutputInfo") << " No EE RecHits.";
606  }
607 
608 
609  //***********************************
610  //*** Get the Jet Rho
611  //***********************************
612  edm::Handle<double> pRho;
613  mEvent.getByLabel(rho_tag_,pRho);
614  if( pRho.isValid() ) {
615  double rho_ = *pRho;
616  if(mRho) mRho->Fill(rho_);
617  }
618 
619  //***********************************
620  //*** Get the Jet collection
621  //***********************************
622  math::XYZTLorentzVector p4tmp[2];
624  mEvent.getByLabel(mInputCollection, pfJets);
625  if (!pfJets.isValid()) return;
626  PFJetCollection::const_iterator jet = pfJets->begin ();
627  int jetIndex = 0;
628  int nJet = 0;
629  int nJetF = 0;
630  int nJetC = 0;
631  for (; jet != pfJets->end (); jet++, jetIndex++) {
632 
633  if (jet->pt() > 10.) {
634  if (fabs(jet->eta()) > 1.5)
635  nJetF++;
636  else
637  nJetC++;
638  }
639  if (jet->pt() > 10.) {
640  if (mEta) mEta->Fill (jet->eta());
641  if (mEtaFineBin) mEtaFineBin->Fill (jet->eta());
642  //if (mEtaFineBin1p) mEtaFineBin1p->Fill (jet->eta());
643  //if (mEtaFineBin2p) mEtaFineBin2p->Fill (jet->eta());
644  //if (mEtaFineBin3p) mEtaFineBin3p->Fill (jet->eta());
645  //if (mEtaFineBin1m) mEtaFineBin1m->Fill (jet->eta());
646  //if (mEtaFineBin2m) mEtaFineBin2m->Fill (jet->eta());
647  //if (mEtaFineBin3m) mEtaFineBin3m->Fill (jet->eta());
648  if (mPhiFineBin) mPhiFineBin->Fill (jet->phi());
649  }
650  if (mjetArea) mjetArea->Fill(jet->jetArea());
651  if (mPhi) mPhi->Fill (jet->phi());
652  if (mE) mE->Fill (jet->energy());
653  if (mE_80) mE_80->Fill (jet->energy());
654  if (mP) mP->Fill (jet->p());
655  if (mP_80) mP_80->Fill (jet->p());
656  if (mPt) mPt->Fill (jet->pt());
657  if (mPt_80) mPt_80->Fill (jet->pt());
658  if (mMass) mMass->Fill (jet->mass());
659  if (mMass_80) mMass_80->Fill (jet->mass());
660  if (mConstituents) mConstituents->Fill (jet->nConstituents());
661  if (mConstituents_80) mConstituents_80->Fill (jet->nConstituents());
662  if (jet == pfJets->begin ()) { // first jet
663  if (mEtaFirst) mEtaFirst->Fill (jet->eta());
664  if (mPhiFirst) mPhiFirst->Fill (jet->phi());
665  if (mPtFirst) mPtFirst->Fill (jet->pt());
666  if (mPtFirst_80) mPtFirst_80->Fill (jet->pt());
667  if (mPtFirst_3000) mPtFirst_3000->Fill (jet->pt());
668  }
669  if (jetIndex == 0) {
670  nJet++;
671  p4tmp[0] = jet->p4();
672  }
673  if (jetIndex == 1) {
674  nJet++;
675  p4tmp[1] = jet->p4();
676  }
677 
678  // if (mMaxEInEmTowers) mMaxEInEmTowers->Fill (jet->maxEInEmTowers());
679  // if (mMaxEInHadTowers) mMaxEInHadTowers->Fill (jet->maxEInHadTowers());
680  // if (mHadEnergyInHO) mHadEnergyInHO->Fill (jet->hadEnergyInHO());
681  // if (mHadEnergyInHO_80) mHadEnergyInHO_80->Fill (jet->hadEnergyInHO());
682  // if (mHadEnergyInHO_3000) mHadEnergyInHO_3000->Fill (jet->hadEnergyInHO());
683  // if (mHadEnergyInHB) mHadEnergyInHB->Fill (jet->hadEnergyInHB());
684  // if (mHadEnergyInHB_80) mHadEnergyInHB_80->Fill (jet->hadEnergyInHB());
685  // if (mHadEnergyInHB_3000) mHadEnergyInHB_3000->Fill (jet->hadEnergyInHB());
686  if (mHadEnergyInHF) mHadEnergyInHF->Fill (jet->HFHadronEnergy());
687  if (mHadEnergyInHF_80) mHadEnergyInHF_80->Fill (jet->HFHadronEnergy());
688  if (mHadEnergyInHF_3000) mHadEnergyInHF_3000->Fill (jet->HFHadronEnergy());
689  // if (mHadEnergyInHE) mHadEnergyInHE->Fill (jet->hadEnergyInHE());
690  // if (mHadEnergyInHE_80) mHadEnergyInHE_80->Fill (jet->hadEnergyInHE());
691  // if (mHadEnergyInHE_3000) mHadEnergyInHE_3000->Fill (jet->hadEnergyInHE());
692  // if (mEmEnergyInEB) mEmEnergyInEB->Fill (jet->emEnergyInEB());
693  // if (mEmEnergyInEB_80) mEmEnergyInEB_80->Fill (jet->emEnergyInEB());
694  // if (mEmEnergyInEB_3000) mEmEnergyInEB_3000->Fill (jet->emEnergyInEB());
695  // if (mEmEnergyInEE) mEmEnergyInEE->Fill (jet->emEnergyInEE());
696  // if (mEmEnergyInEE_80) mEmEnergyInEE_80->Fill (jet->emEnergyInEE());
697  // if (mEmEnergyInEE_3000) mEmEnergyInEE_3000->Fill (jet->emEnergyInEE());
698  if (mEmEnergyInHF) mEmEnergyInHF->Fill (jet->HFEMEnergy());
699  if (mEmEnergyInHF_80) mEmEnergyInHF_80->Fill (jet->HFEMEnergy());
700  if (mEmEnergyInHF_3000) mEmEnergyInHF_3000->Fill (jet->HFEMEnergy());
701  // if (mEnergyFractionHadronic) mEnergyFractionHadronic->Fill (jet->energyFractionHadronic());
702  // if (mEnergyFractionEm) mEnergyFractionEm->Fill (jet->emEnergyFraction());
703  //
704  // if (mHFTotal) mHFTotal->Fill (jet->hadEnergyInHF()+jet->emEnergyInHF());
705  // if (mHFTotal_80) mHFTotal_80->Fill (jet->hadEnergyInHF()+jet->emEnergyInHF());
706  // if (mHFTotal_3000) mHFTotal_3000->Fill (jet->hadEnergyInHF()+jet->emEnergyInHF());
707  // if (mHFLong) mHFLong->Fill (jet->hadEnergyInHF()*0.5+jet->emEnergyInHF());
708  // if (mHFLong_80) mHFLong_80->Fill (jet->hadEnergyInHF()*0.5+jet->emEnergyInHF());
709  // if (mHFLong_3000) mHFLong_3000->Fill (jet->hadEnergyInHF()*0.5+jet->emEnergyInHF());
710  // if (mHFShort) mHFShort->Fill (jet->hadEnergyInHF()*0.5);
711  // if (mHFShort_80) mHFShort_80->Fill (jet->hadEnergyInHF()*0.5);
712  // if (mHFShort_3000) mHFShort_3000->Fill (jet->hadEnergyInHF()*0.5);
713 
714  if (mChargedEmEnergy_80) mChargedEmEnergy_80->Fill (jet->chargedEmEnergy());
715  if (mChargedHadronEnergy_80) mChargedHadronEnergy_80->Fill (jet->chargedHadronEnergy());
716  if (mNeutralEmEnergy_80) mNeutralEmEnergy_80->Fill (jet->neutralEmEnergy());
717  if (mNeutralHadronEnergy_80) mNeutralHadronEnergy_80->Fill (jet->neutralHadronEnergy());
718 
719  if (mChargedEmEnergy_3000) mChargedEmEnergy_3000->Fill (jet->chargedEmEnergy());
720  if (mChargedHadronEnergy_3000) mChargedHadronEnergy_3000->Fill (jet->chargedHadronEnergy());
721  if (mNeutralEmEnergy_3000) mNeutralEmEnergy_3000->Fill (jet->neutralEmEnergy());
722  if (mNeutralHadronEnergy_3000) mNeutralHadronEnergy_3000->Fill (jet->neutralHadronEnergy());
723 
724  if (fabs(jet->eta())<1.5) mChargedEmEnergyFraction_B->Fill (jet->chargedEmEnergyFraction());
725  if (fabs(jet->eta())>1.5 && fabs(jet->eta())<3.0) mChargedEmEnergyFraction_E->Fill (jet->chargedEmEnergyFraction());
726  if (fabs(jet->eta())>3.0 && fabs(jet->eta())<6.0) mChargedEmEnergyFraction_F->Fill (jet->chargedEmEnergyFraction());
727  if (fabs(jet->eta())<1.5) mChargedHadronEnergyFraction_B->Fill (jet->chargedHadronEnergyFraction());
728  if (fabs(jet->eta())>1.5 && fabs(jet->eta())<3.0) mChargedHadronEnergyFraction_E->Fill (jet->chargedHadronEnergyFraction());
729  if (fabs(jet->eta())>3.0 && fabs(jet->eta())<6.0) mChargedHadronEnergyFraction_F->Fill (jet->chargedHadronEnergyFraction());
730  if (fabs(jet->eta())<1.5) mNeutralEmEnergyFraction_B->Fill (jet->neutralEmEnergyFraction());
731  if (fabs(jet->eta())>1.5 && fabs(jet->eta())<3.0) mNeutralEmEnergyFraction_E->Fill (jet->neutralEmEnergyFraction());
732  if (fabs(jet->eta())>3.0 && fabs(jet->eta())<6.0) mNeutralEmEnergyFraction_F->Fill (jet->neutralEmEnergyFraction());
733  if (fabs(jet->eta())<1.5) mNeutralHadronEnergyFraction_B->Fill (jet->neutralHadronEnergyFraction());
734  if (fabs(jet->eta())>1.5 && fabs(jet->eta())<3.0) mNeutralHadronEnergyFraction_E->Fill (jet->neutralHadronEnergyFraction());
735  if (fabs(jet->eta())>3.0 && fabs(jet->eta())<6.0) mNeutralHadronEnergyFraction_F->Fill (jet->neutralHadronEnergyFraction());
736 
737  // if (mN90) mN90->Fill (jet->n90());
738  //mJetEnergyProfile->Fill (jet->eta(), jet->phi(), jet->energy());
739  // mHadJetEnergyProfile->Fill (jet->eta(), jet->phi(), jet->hadEnergyInHO()+jet->hadEnergyInHB()+jet->hadEnergyInHF()+jet->hadEnergyInHE());
740  // mEMJetEnergyProfile->Fill (jet->eta(), jet->phi(), jet->emEnergyInEB()+jet->emEnergyInEE()+jet->emEnergyInHF());
741  }
742 
743 
744 
745 
746  if (mNJetsEtaC) mNJetsEtaC->Fill( nJetC );
747  if (mNJetsEtaF) mNJetsEtaF->Fill( nJetF );
748 
749  if (nJet == 2) {
750  if (mMjj) mMjj->Fill( (p4tmp[0]+p4tmp[1]).mass() );
751  if (mMjj_3000) mMjj_3000->Fill( (p4tmp[0]+p4tmp[1]).mass() );
752  }
753 
754  // Count Jets above Pt cut
755  for (int istep = 0; istep < 100; ++istep) {
756  int njet = 0;
757  float ptStep = (istep * (200./100.));
758 
759  for ( PFJetCollection::const_iterator cal = pfJets->begin(); cal != pfJets->end(); ++ cal ) {
760  if ( cal->pt() > ptStep ) njet++;
761  }
762  mNJets1->Fill( ptStep, njet );
763  }
764 
765  for (int istep = 0; istep < 100; ++istep) {
766  int njet = 0;
767  float ptStep = (istep * (4000./100.));
768  for ( PFJetCollection::const_iterator cal = pfJets->begin(); cal != pfJets->end(); ++ cal ) {
769  if ( cal->pt() > ptStep ) njet++;
770  }
771  mNJets2->Fill( ptStep, njet );
772  }
773 
774 if (!mEvent.isRealData()){
775  // Gen jet analysis
776  Handle<GenJetCollection> genJets;
777  mEvent.getByLabel(mInputGenCollection, genJets);
778  if (!genJets.isValid()) return;
779  GenJetCollection::const_iterator gjet = genJets->begin ();
780  int gjetIndex = 0;
781  for (; gjet != genJets->end (); gjet++, gjetIndex++) {
782  if (mGenEta) mGenEta->Fill (gjet->eta());
783  if (mGenPhi) mGenPhi->Fill (gjet->phi());
784  if (mGenPt) mGenPt->Fill (gjet->pt());
785  if (mGenPt_80) mGenPt_80->Fill (gjet->pt());
786  if (gjet == genJets->begin ()) { // first jet
787  if (mGenEtaFirst) mGenEtaFirst->Fill (gjet->eta());
788  if (mGenPhiFirst) mGenPhiFirst->Fill (gjet->phi());
789  }
790  }
791 
792 
793  // now match PFJets to GenJets
794  JetMatchingTools jetMatching (mEvent);
795  if (!(mInputGenCollection.label().empty())) {
796  // Handle<GenJetCollection> genJets;
797  // mEvent.getByLabel(mInputGenCollection, genJets);
798 
799  std::vector <std::vector <const reco::GenParticle*> > genJetConstituents (genJets->size());
800  std::vector <std::vector <const reco::GenParticle*> > pfJetConstituents (pfJets->size());
801  /*
802  if (mRThreshold > 0) {
803  }
804  else {
805  for (unsigned iGenJet = 0; iGenJet < genJets->size(); ++iGenJet) {
806  genJetConstituents [iGenJet] = jetMatching.getGenParticles ((*genJets) [iGenJet]);
807  }
808 
809  for (unsigned iPFJet = 0; iPFJet < pfJets->size(); ++iPFJet) {
810  pfJetConstituents [iPFJet] = jetMatching.getGenParticles ((*pfJets) [iPFJet], false);
811  }
812  }
813  */
814 
815  for (unsigned iGenJet = 0; iGenJet < genJets->size(); ++iGenJet) { //****************************************************************
816  //for (unsigned iGenJet = 0; iGenJet < 1; ++iGenJet) { // only FIRST Jet !!!!
817  const GenJet& genJet = (*genJets) [iGenJet];
818  double genJetPt = genJet.pt();
819 
820  //std::cout << iGenJet <<". Genjet: pT = " << genJetPt << "GeV" << std::endl; // *****************************************************
821 
822  if (fabs(genJet.eta()) > 6.) continue; // out of detector
823  if (genJetPt < mMatchGenPtThreshold) continue; // no low momentum
824  //double logPtGen = log10 (genJetPt);
825  //mAllGenJetsPt->Fill (logPtGen);
826  //mAllGenJetsEta->Fill (logPtGen, genJet.eta());
827  if (pfJets->size() <= 0) continue; // no PFJets - nothing to match
828  if (mRThreshold > 0) {
829  unsigned iPFJetBest = 0;
830  double deltaRBest = 999.;
831  for (unsigned iPFJet = 0; iPFJet < pfJets->size(); ++iPFJet) {
832  double dR = deltaR (genJet.eta(), genJet.phi(), (*pfJets) [iPFJet].eta(), (*pfJets) [iPFJet].phi());
833  if (deltaRBest < mRThreshold && dR < mRThreshold && genJet.pt() > 5.) {
834  /*
835  std::cout << "Yet another matched jet for GenJet pt=" << genJet.pt()
836  << " previous PFJet pt/dr: " << (*pfJets) [iPFJetBest].pt() << '/' << deltaRBest
837  << " new PFJet pt/dr: " << (*pfJets) [iPFJet].pt() << '/' << dR
838  << std::endl;
839  */
840  }
841  if (dR < deltaRBest) {
842  iPFJetBest = iPFJet;
843  deltaRBest = dR;
844  }
845  }
846  if (mTurnOnEverything.compare("yes")==0) {
847  //mRMatch->Fill (logPtGen, genJet.eta(), deltaRBest);
848  }
849  if (deltaRBest < mRThreshold) { // Matched
850  fillMatchHists (genJet, (*pfJets) [iPFJetBest]);
851  }
852  }
853  /*
854  else {
855  unsigned iPFJetBest = 0;
856  double energyFractionBest = 0.;
857  for (unsigned iPFJet = 0; iPFJet < pfJets->size(); ++iPFJet) {
858  double energyFraction = jetMatching.overlapEnergyFraction (genJetConstituents [iGenJet],
859  pfJetConstituents [iPFJet]);
860  if (energyFraction > energyFractionBest) {
861  iPFJetBest = iPFJet;
862  energyFractionBest = energyFraction;
863  }
864  }
865  if (mTurnOnEverything.compare("yes")==0) {
866  mGenJetMatchEnergyFraction->Fill (logPtGen, genJet.eta(), energyFractionBest);
867  }
868  if (energyFractionBest > mGenEnergyFractionThreshold) { // good enough
869  double reverseEnergyFraction = jetMatching.overlapEnergyFraction (pfJetConstituents [iPFJetBest],
870  genJetConstituents [iGenJet]);
871  if (mTurnOnEverything.compare("yes")==0) {
872  mReverseMatchEnergyFraction->Fill (logPtGen, genJet.eta(), reverseEnergyFraction);
873  }
874  if (reverseEnergyFraction > mReverseEnergyFractionThreshold) { // Matched
875  fillMatchHists (genJet, (*pfJets) [iPFJetBest]);
876  }
877  }
878  }
879  */
880  }
881  }
882 }
883 
884 }
885 
886 void PFJetTesterUnCorr::fillMatchHists (const reco::GenJet& fGenJet, const reco::PFJet& fPFJet) {
887  double logPtGen = log10 (fGenJet.pt());
888  double PtGen = fGenJet.pt();
889  double PtPF = fPFJet.pt();
890  //mMatchedGenJetsPt->Fill (logPtGen);
891  //mMatchedGenJetsEta->Fill (logPtGen, fGenJet.eta());
892 
893  double PtThreshold = 10.;
894 
895  if (mTurnOnEverything.compare("yes")==0) {
896  mDeltaEta->Fill (logPtGen, fGenJet.eta(), fPFJet.eta()-fGenJet.eta());
897  mDeltaPhi->Fill (logPtGen, fGenJet.eta(), fPFJet.phi()-fGenJet.phi());
898  //mEScale->Fill (logPtGen, fGenJet.eta(), fPFJet.energy()/fGenJet.energy());
899  //mlinEScale->Fill (fGenJet.pt(), fGenJet.eta(), fPFJet.energy()/fGenJet.energy());
900  //mDeltaE->Fill (logPtGen, fGenJet.eta(), fPFJet.energy()-fGenJet.energy());
901 
902  mEScaleFineBin->Fill (logPtGen, fGenJet.eta(), fPFJet.energy()/fGenJet.energy());
903 
904  if (fGenJet.pt()>PtThreshold) {
905  mEScale_pt10->Fill (logPtGen, fGenJet.eta(), fPFJet.energy()/fGenJet.energy());
906 
907  }
908 
909  }
910  if (fPFJet.pt() > PtThreshold) {
911  mDelEta->Fill (fGenJet.eta()-fPFJet.eta());
912  mDelPhi->Fill (fGenJet.phi()-fPFJet.phi());
913  mDelPt->Fill ((fGenJet.pt()-fPFJet.pt())/fGenJet.pt());
914  }
915 
916  if (fabs(fGenJet.eta())<1.5) {
917 
918  //mpTScaleB_s->Fill (log10(PtGen), PtPF/PtGen);
919  mpTScaleB_d->Fill (log10(PtGen), PtPF/PtGen);
920  mpTScalePhiB_d->Fill (fGenJet.phi(), PtPF/PtGen);
921 
922  if (PtGen>30.0 && PtGen<200.0) {
923  mpTScale1DB_30_200->Fill (fPFJet.pt()/fGenJet.pt());
924  }
925  if (PtGen>200.0 && PtGen<600.0) {
926  mpTScale1DB_200_600->Fill (fPFJet.pt()/fGenJet.pt());
927  }
928  if (PtGen>600.0 && PtGen<1500.0) {
929  mpTScale1DB_600_1500->Fill (fPFJet.pt()/fGenJet.pt());
930  }
931  if (PtGen>1500.0 && PtGen<3500.0) {
932  mpTScale1DB_1500_3500->Fill (fPFJet.pt()/fGenJet.pt());
933  }
934 
935  }
936 
937  if (fabs(fGenJet.eta())>1.5 && fabs(fGenJet.eta())<3.0) {
938 
939  //mpTScaleE_s->Fill (log10(PtGen), PtPF/PtGen);
940  mpTScaleE_d->Fill (log10(PtGen), PtPF/PtGen);
941  mpTScalePhiE_d->Fill (fGenJet.phi(), PtPF/PtGen);
942 
943  if (PtGen>30.0 && PtGen<200.0) {
944  mpTScale1DE_30_200->Fill (fPFJet.pt()/fGenJet.pt());
945  }
946  if (PtGen>200.0 && PtGen<600.0) {
947  mpTScale1DE_200_600->Fill (fPFJet.pt()/fGenJet.pt());
948  }
949  if (PtGen>600.0 && PtGen<1500.0) {
950  mpTScale1DE_600_1500->Fill (fPFJet.pt()/fGenJet.pt());
951  }
952  if (PtGen>1500.0 && PtGen<3500.0) {
953  mpTScale1DE_1500_3500->Fill (fPFJet.pt()/fGenJet.pt());
954  }
955 
956  }
957 
958  if (fabs(fGenJet.eta())>3.0 && fabs(fGenJet.eta())<6.0) {
959 
960  //mpTScaleF_s->Fill (log10(PtGen), PtPF/PtGen);
961  mpTScaleF_d->Fill (log10(PtGen), PtPF/PtGen);
962  mpTScalePhiF_d->Fill (fGenJet.phi(), PtPF/PtGen);
963 
964  if (PtGen>30.0 && PtGen<200.0) {
965  mpTScale1DF_30_200->Fill (fPFJet.pt()/fGenJet.pt());
966  }
967  if (PtGen>200.0 && PtGen<600.0) {
968  mpTScale1DF_200_600->Fill (fPFJet.pt()/fGenJet.pt());
969  }
970  if (PtGen>600.0 && PtGen<1500.0) {
971  mpTScale1DF_600_1500->Fill (fPFJet.pt()/fGenJet.pt());
972  }
973  if (PtGen>1500.0 && PtGen<3500.0) {
974  mpTScale1DF_1500_3500->Fill (fPFJet.pt()/fGenJet.pt());
975  }
976 
977  }
978 
979  if (fGenJet.pt()>30.0 && fGenJet.pt()<200.0) {
980  //mpTScale_30_200_s->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
981  mpTScale_30_200_d->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
982  //mpTScale1D_30_200->Fill (fPFJet.pt()/fGenJet.pt());
983  }
984 
985  if (fGenJet.pt()>200.0 && fGenJet.pt()<600.0) {
986  //mpTScale_200_600_s->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
987  mpTScale_200_600_d->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
988  //mpTScale1D_200_600->Fill (fPFJet.pt()/fGenJet.pt());
989  }
990 
991  if (fGenJet.pt()>600.0 && fGenJet.pt()<1500.0) {
992  //mpTScale_600_1500_s->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
993  mpTScale_600_1500_d->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
994  //mpTScale1D_600_1500->Fill (fPFJet.pt()/fGenJet.pt());
995  }
996 
997  if (fGenJet.pt()>1500.0 && fGenJet.pt()<3500.0) {
998  //mpTScale_1500_3500_s->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
999  mpTScale_1500_3500_d->Fill (fGenJet.eta(),fPFJet.pt()/fGenJet.pt());
1000  //mpTScale1D_1500_3500->Fill (fPFJet.pt()/fGenJet.pt());
1001  }
1002 
1003 
1004 
1005 }
1006 
MonitorElement * mDeltaPhi
MonitorElement * mPhiFineBin
MonitorElement * mChargedHadronEnergyFraction_F
void getManyByType(std::vector< Handle< PROD > > &results) const
Definition: Event.h:408
MonitorElement * mNJetsEtaC
int i
Definition: DBlmapReader.cc:9
MonitorElement * mNeutralEmEnergyFraction_F
MonitorElement * mEmEnergyInHF_3000
MonitorElement * mHadTiming
MonitorElement * numberofevents
MonitorElement * mGenPhiFirst
MonitorElement * mpTScale1DE_1500_3500
MonitorElement * mNeutralEmEnergyFraction_B
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:717
MonitorElement * mMass_80
MonitorElement * mNJetsEtaF
MonitorElement * mNeutralHadronEnergyFraction_E
MonitorElement * mpTScalePhiE_d
Base class for all types of Jets.
Definition: Jet.h:21
MonitorElement * mMjj_3000
std::vector< CaloTower >::const_iterator const_iterator
MonitorElement * mDelEta
MonitorElement * mPtFirst
MonitorElement * mEtaFirst
MonitorElement * mpTScale1DE_600_1500
MonitorElement * mChargedEmEnergy_3000
MonitorElement * mpTScaleB_d
MonitorElement * mGenEtaFirst
PFJetTesterUnCorr(const edm::ParameterSet &)
MonitorElement * mNeutralEmEnergyFraction_E
std::string mOutputFile
int njet
Definition: HydjetWrapper.h:91
bool isRealData() const
Definition: EventBase.h:60
Jets made from PFObjects.
Definition: PFJet.h:22
MonitorElement * mP
MonitorElement * mE_80
virtual double eta() const
momentum pseudorapidity
edm::InputTag mInputGenCollection
MonitorElement * mChargedHadronEnergyFraction_B
MonitorElement * mChargedHadronEnergyFraction_E
void Fill(long long x)
virtual void endJob()
MonitorElement * mNeutralHadronEnergy_80
MonitorElement * mEta
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
MonitorElement * mDelPhi
MonitorElement * mpTScale1DF_200_600
MonitorElement * mPthat_3000
virtual double energy() const
energy
MonitorElement * mGenPt_80
edm::InputTag rho_tag_
MonitorElement * mEmTiming
MonitorElement * mpTScale1DB_600_1500
MonitorElement * mpTScale_1500_3500_d
MonitorElement * mpTScale1DE_200_600
MonitorElement * mpTScaleE_d
MonitorElement * mPt_80
MonitorElement * mNeutralHadronEnergy_3000
MonitorElement * mpTScale1DF_1500_3500
MonitorElement * mP_80
std::string mTurnOnEverything
edm::InputTag mInputCollection
MonitorElement * mpTScale_600_1500_d
MonitorElement * mPhiFirst
int j
Definition: DBlmapReader.cc:9
void fillMatchHists(const reco::GenJet &fGenJet, const reco::PFJet &fPFJet)
Gen Close.
Jets made from MC generator particles.
Definition: GenJet.h:25
MonitorElement * mChargedHadronEnergy_80
MonitorElement * bookProfile(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, const char *option="s")
Definition: DQMStore.cc:1031
MonitorElement * mChargedEmEnergyFraction_E
MonitorElement * mMass
MonitorElement * mpTScalePhiF_d
MonitorElement * mGenPt
bool isValid() const
Definition: HandleBase.h:76
MonitorElement * mEmEnergyInHF
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
MonitorElement * mPthat_80
MonitorElement * mpTScale1DB_30_200
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
MonitorElement * mNeutralEmEnergy_80
MonitorElement * mE
MonitorElement * mRho
MonitorElement * mpTScaleF_d
MonitorElement * mDeltaEta
MonitorElement * mEScale_pt10
MonitorElement * mNeutralHadronEnergyFraction_F
MonitorElement * mPtFirst_3000
MonitorElement * mChargedHadronEnergy_3000
MonitorElement * mEtaFineBin
virtual double pt() const
transverse momentum
virtual void beginJob()
MonitorElement * mPhi
MonitorElement * mpTScale1DB_200_600
MonitorElement * mConstituents
MonitorElement * mChargedEmEnergyFraction_B
std::string const & label() const
Definition: InputTag.h:25
MonitorElement * mNeutralHadronEnergyFraction_B
MonitorElement * mpTScalePhiB_d
MonitorElement * mpTScale_30_200_d
MonitorElement * mpTScale1DF_30_200
MonitorElement * mjetArea
MonitorElement * mGenPhi
MonitorElement * mConstituents_80
MonitorElement * mChargedEmEnergy_80
tuple mass
Definition: scaleCards.py:27
MonitorElement * mpTScale1DB_1500_3500
MonitorElement * mGenEta
MonitorElement * mpTScale_200_600_d
MonitorElement * mMjj
MonitorElement * mHadEnergyInHF_80
MonitorElement * mHadEnergyInHF
MonitorElement * mNeutralEmEnergy_3000
MonitorElement * mChargedEmEnergyFraction_F
tuple pfJets
Definition: pfJets_cff.py:8
virtual double phi() const
momentum azimuthal angle
MonitorElement * mPtFirst_80
MonitorElement * mpTScale1DF_600_1500
MonitorElement * mNJets2
MonitorElement * mEScaleFineBin
MonitorElement * mpTScale1DE_30_200
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:429
MonitorElement * mHadEnergyInHF_3000
MonitorElement * mNJets1
virtual void analyze(const edm::Event &, const edm::EventSetup &)
MonitorElement * mEmEnergyInHF_80
MonitorElement * mDelPt
MonitorElement * mPt