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

Comparing UserCode/CmsHi/JetAnalysis/src/HiJetResponseAnalyzer.cc (file contents):
Revision 1.3 by yilmaz, Fri Sep 10 14:47:48 2010 UTC vs.
Revision 1.6 by yilmaz, Mon Oct 18 16:13:37 2010 UTC

# Line 42 | Line 42
42  
43   #include "CommonTools/UtilAlgos/interface/TFileService.h"
44   #include "FWCore/ServiceRegistry/interface/Service.h"
45 <
45 > #include "DataFormats/Math/interface/deltaR.h"
46  
47   #include "TTree.h"
48  
# Line 62 | Line 62 | struct JRA{
62     float b;
63     float hf;
64     float jtpt[MAXJETS];
65 +   float jtcorpt[MAXJETS];
66     float refpt[MAXJETS];
67     float jteta[MAXJETS];
68     float refeta[MAXJETS];
# Line 71 | Line 72 | struct JRA{
72     float weight;
73   };
74  
75 +
76 + struct JRAV{
77 +
78 +  float jtpt;
79 +  float jtcorpt;
80 +  float refpt;
81 +  float refcorpt;
82 +  float jteta;
83 +  float refeta;
84 +  float jtphi;
85 +  float refphi;
86 +
87 + };
88 +
89   //
90   // class declaration
91   //
92 + bool compareCorPt(JRAV a, JRAV b) {return a.jtcorpt > b.jtcorpt;}
93 + bool comparePt(JRAV a, JRAV b) {return a.jtpt > b.jtpt;}
94  
95   class HiJetResponseAnalyzer : public edm::EDAnalyzer {
96     public:
# Line 85 | Line 102 | class HiJetResponseAnalyzer : public edm
102        virtual void beginJob() ;
103        virtual void analyze(const edm::Event&, const edm::EventSetup&);
104        virtual void endJob() ;
105 <
105 >   bool selectJet(int i);
106        // ----------member data ---------------------------
107  
108 <   bool usePat_;
109 <   bool doMC_;
110 <
111 <   double genJetPtMin_;
112 <
108 >  bool usePat_;
109 >  bool doMC_;
110 >  bool filterJets_;
111 >  bool diJetsOnly_;
112 >  bool matchDiJets_;
113 >  bool matchPatGen_;
114 >  bool matchNew_;
115 >  bool sortJets_;
116 >  bool correctJets_;
117 >
118 >  double matchR_;  
119 >   double genPtMin_;
120 >   double ptMin_;
121 >   double emfMin_;
122 >   double n90Min_;
123 >   double n90hitMin_;
124  
125     edm::InputTag jetTag_;
126 +  edm::InputTag matchTag_;
127  
128     JRA jra_;
129     TTree* t;
# Line 103 | Line 132 | class HiJetResponseAnalyzer : public edm
132     edm::Handle<reco::Centrality> cent;
133  
134     edm::Handle<reco::JetView> jets;
135 <   //   edm::Handle<pat::JetCollection> jets;
136 <
135 >   edm::Handle<pat::JetCollection> patjets;
136 >  edm::Handle<reco::JetView> matchedJets;
137  
138     edm::Service<TFileService> fs;
139  
140   };
141  
142 + bool HiJetResponseAnalyzer::selectJet(int i){
143 +   const reco::Jet& jet = (*jets)[i];
144 +   if(usePat_){
145 +      const pat::Jet& patjet = (*patjets)[i];
146 +      if(patjet.emEnergyFraction() <= emfMin_) return false;
147 +      if(patjet.jetID().n90Hits <= n90hitMin_) return false;
148 +      if(doMC_){
149 +        
150 +      }
151 +
152 +   }
153 +
154 +   return true;
155 + }
156 +
157 +
158   //
159   // constants, enums and typedefs
160   //
# Line 125 | Line 170 | HiJetResponseAnalyzer::HiJetResponseAnal
170  
171   {
172     //now do what ever initialization is needed
173 <   genJetPtMin_ = 0;
173 >  matchR_ = iConfig.getUntrackedParameter<double>("matchR",0.25);
174 >
175 >   ptMin_ = iConfig.getUntrackedParameter<double>("ptMin",0);
176 >   genPtMin_ = iConfig.getUntrackedParameter<double>("genPtMin",20);
177 >   emfMin_ = iConfig.getUntrackedParameter<double>("emfMin",0.01);
178 >   n90Min_ = iConfig.getUntrackedParameter<double>("n90Min",1);
179 >   n90hitMin_ = iConfig.getUntrackedParameter<double>("n90hitMin",1);
180  
181 <   usePat_ = false;
182 <   doMC_ = false;
181 >   filterJets_ = iConfig.getUntrackedParameter<bool>("filterJets",true);
182 >   diJetsOnly_ = iConfig.getUntrackedParameter<bool>("diJetsOnly",false);
183 >   matchDiJets_ = iConfig.getUntrackedParameter<bool>("matchDiJets",false);
184 >   matchPatGen_ = iConfig.getUntrackedParameter<bool>("matchPatGen",false);
185  
186 <   jetTag_ = edm::InputTag("selectedPatJets");
186 >   matchNew_ = iConfig.getUntrackedParameter<bool>("matchNew",false);
187  
188 +   usePat_ = iConfig.getUntrackedParameter<bool>("usePat",true);
189 +   doMC_ = iConfig.getUntrackedParameter<bool>("doMC",true);
190 +
191 +   sortJets_ = iConfig.getUntrackedParameter<bool>("sortJets",true);
192 +   correctJets_ = iConfig.getUntrackedParameter<bool>("correctJets",false);
193 +
194 +   jetTag_ = iConfig.getUntrackedParameter<edm::InputTag>("src",edm::InputTag("selectedPatJets"));
195 +   matchTag_ = iConfig.getUntrackedParameter<edm::InputTag>("match",edm::InputTag("selectedPatJets"));
196   }
197  
198  
# Line 155 | Line 216 | HiJetResponseAnalyzer::analyze(const edm
216     using namespace edm;
217  
218     iEvent.getByLabel(jetTag_,jets);
219 +   if(usePat_)iEvent.getByLabel(jetTag_,patjets);
220 +   if(matchNew_)iEvent.getByLabel(matchTag_,matchedJets);
221 +
222 +   std::vector<JRAV> jraV;
223  
159   jra_.nref = 0;
224     for(unsigned int j = 0 ; j < jets->size(); ++j){
225 <      
226 <      //      const pat::Jet& jet = (*jets)[j];
227 <      const reco::Jet jet = (*jets)[j];      
228 <      jra_.jtpt[jra_.nref] = jet.pt();
229 <      jra_.jteta[jra_.nref] = jet.eta();
230 <      jra_.jtphi[jra_.nref] = jet.phi();
231 <      
232 <      if(usePat_){
233 <         const pat::Jet& patjet = (*jets)[j];
234 <        
235 <         if(doMC_ && patjet.genJet() != 0){        
236 <            //      if(jet.genJet()->pt() < genJetPtMin_) continue;
237 <            jra_.refpt[jra_.nref] = patjet.genJet()->pt();
238 <            jra_.refeta[jra_.nref] = patjet.genJet()->eta();
239 <            jra_.refphi[jra_.nref] = patjet.genJet()->phi();
240 <                            
241 <         }else{
242 <            jra_.refpt[jra_.nref] = -99;
243 <            jra_.refeta[jra_.nref] = -99;
244 <            jra_.refphi[jra_.nref] = -99;
225 >
226 >     if(filterJets_ && !selectJet(j)) continue;
227 >     const reco::Jet& jet = (*jets)[j];
228 >     JRAV jv;
229 >     jv.jtpt = jet.pt();
230 >     jv.jteta = jet.eta();
231 >     jv.jtphi = jet.phi();
232 >     jv.jtcorpt = jet.pt();
233 >
234 >     if(usePat_){
235 >       const pat::Jet& patjet = (*patjets)[j];
236 >
237 >       jv.jtpt = patjet.correctedJet("raw").pt();
238 >       jv.jtcorpt = patjet.pt();
239 >
240 >       if(doMC_ && matchPatGen_ && patjet.genJet() != 0){
241 >         if(patjet.genJet()->pt() < genPtMin_) continue;
242 >         jv.refpt = patjet.genJet()->pt();
243 >         jv.refeta = patjet.genJet()->eta();
244 >         jv.refphi = patjet.genJet()->phi();
245 >
246 >       }else{
247 >         jv.refpt = -99;
248 >         jv.refeta = -99;
249 >         jv.refphi = -99;
250 >       }
251 >     }
252 >
253 >     if(matchNew_){
254 >       for(unsigned int m = 0 ; m < matchedJets->size(); ++m){
255 >         const reco::Jet& match = (*matchedJets)[m];
256 >         double dr = reco::deltaR(jet.eta(),jet.phi(),match.eta(),match.phi());
257 >         if(dr < matchR_){
258 >           jv.refcorpt = -99;
259 >           jv.refpt = match.pt();
260 >           jv.refeta = match.eta();
261 >           jv.refphi = match.phi();
262           }
263 <      }
264 <      
265 <      jra_.nref++;
266 <      
263 >       }
264 >     }
265 >     jraV.push_back(jv);
266 >   }
267 >
268 >   if(sortJets_){
269 >     if(usePat_ || correctJets_) std::sort(jraV.begin(),jraV.end(),compareCorPt);
270 >     else std::sort(jraV.begin(),jraV.end(),comparePt);
271 >   }
272 >
273 >   for(unsigned int i = 0; i < jraV.size(); ++i){
274 >     JRAV& jv = jraV[i];
275 >     jra_.jtpt[jra_.nref] = jv.jtpt;
276 >     jra_.jteta[jra_.nref] = jv.jteta;
277 >     jra_.jtphi[jra_.nref] = jv.jtphi;
278 >     jra_.jtcorpt[jra_.nref] = jv.jtcorpt;
279 >     jra_.refpt[jra_.nref] = jv.refpt;
280 >     jra_.refeta[jra_.nref] = jv.refeta;
281 >     jra_.refphi[jra_.nref] = jv.refphi;
282     }
283 +   jra_.nref = jraV.size();
284  
285     t->Fill();
286  
# Line 200 | Line 297 | HiJetResponseAnalyzer::beginJob()
297     t->Branch("hf",&jra_.hf,"hf/F");
298     t->Branch("nref",&jra_.nref,"nref/I");
299     t->Branch("jtpt",jra_.jtpt,"jtpt[nref]/F");
300 +   t->Branch("jtcorpt",jra_.jtcorpt,"jtcorpt[nref]/F");
301     t->Branch("refpt",jra_.refpt,"refpt[nref]/F");
302 +   t->Branch("refcorpt",jra_.refpt,"refcorpt[nref]/F");
303     t->Branch("jteta",jra_.jteta,"jteta[nref]/F");
304     t->Branch("refeta",jra_.refeta,"refeta[nref]/F");
305     t->Branch("jtphi",jra_.jtphi,"jtphi[nref]/F");

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines