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.5 by yilmaz, Wed Sep 29 15:59:10 2010 UTC vs.
Revision 1.9 by yilmaz, Sun Oct 24 14:39:32 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 72 | Line 72 | struct JRA{
72     float weight;
73   };
74  
75 + struct JRAV{
76 +  int index;
77 +  float jtpt;
78 +  float jtcorpt;
79 +  float refpt;
80 +  float refcorpt;
81 +  float jteta;
82 +  float refeta;
83 +  float jtphi;
84 +  float refphi;
85 +
86 + };
87 +
88   //
89   // class declaration
90   //
91 + bool compareCorPt(JRAV a, JRAV b) {return a.jtcorpt > b.jtcorpt;}
92 + bool comparePt(JRAV a, JRAV b) {return a.jtpt > b.jtpt;}
93  
94   class HiJetResponseAnalyzer : public edm::EDAnalyzer {
95     public:
# Line 87 | Line 102 | class HiJetResponseAnalyzer : public edm
102        virtual void analyze(const edm::Event&, const edm::EventSetup&);
103        virtual void endJob() ;
104     bool selectJet(int i);
90 –
105        // ----------member data ---------------------------
106  
107 <   bool usePat_;
108 <   bool doMC_;
109 <   bool filterJets_;
110 <   bool diJetsOnly_;
111 <   bool matchDiJets_;
107 >  bool usePat_;
108 >  bool doMC_;
109 >  bool filterJets_;
110 >  bool diJetsOnly_;
111 >  bool matchDiJets_;
112 >  bool matchPatGen_;
113 >  bool matchNew_;
114 >  bool sortJets_;
115 >  bool correctJets_;
116  
117 +  double matchR_;  
118     double genPtMin_;
119     double ptMin_;
120     double emfMin_;
121     double n90Min_;
122     double n90hitMin_;
123  
124 <   edm::InputTag jetTag_;
124 >  edm::InputTag jetTag_;
125 >  edm::InputTag matchTag_;
126 >  std::vector<edm::InputTag> matchTags_;
127 >  
128 >  JRA jra_;
129 >  std::vector<JRA> jraMatch_;
130  
107 –   JRA jra_;
131     TTree* t;
132  
133     edm::Handle<edm::GenHIEvent> mc;
# Line 112 | Line 135 | class HiJetResponseAnalyzer : public edm
135  
136     edm::Handle<reco::JetView> jets;
137     edm::Handle<pat::JetCollection> patjets;
138 +  edm::Handle<reco::JetView> matchedJets;
139  
140     edm::Service<TFileService> fs;
141  
142   };
143  
144   bool HiJetResponseAnalyzer::selectJet(int i){
121 –  
145     const reco::Jet& jet = (*jets)[i];
123 –  
146     if(usePat_){
125 –      //      cout<<"a"<<endl;
147        const pat::Jet& patjet = (*patjets)[i];
127 –      //      cout<<"b"<<endl;
128 –      
148        if(patjet.emEnergyFraction() <= emfMin_) return false;
130 –      //      cout<<"c"<<endl;
131 –
149        if(patjet.jetID().n90Hits <= n90hitMin_) return false;
133 –      //      cout<<"d"<<endl;
134 –
150        if(doMC_){
151          
152        }
# Line 139 | Line 154 | bool HiJetResponseAnalyzer::selectJet(in
154     }
155  
156     return true;
142 –
157   }
158  
159  
146 –
160   //
161   // constants, enums and typedefs
162   //
# Line 159 | Line 172 | HiJetResponseAnalyzer::HiJetResponseAnal
172  
173   {
174     //now do what ever initialization is needed
175 +  matchR_ = iConfig.getUntrackedParameter<double>("matchR",0.25);
176  
177     ptMin_ = iConfig.getUntrackedParameter<double>("ptMin",0);
178 <   genPtMin_ = iConfig.getUntrackedParameter<double>("genPtMin",0);
179 <   emfMin_ = iConfig.getUntrackedParameter<double>("emfMin",0);
180 <   n90Min_ = iConfig.getUntrackedParameter<double>("n90Min",0);
181 <   n90hitMin_ = iConfig.getUntrackedParameter<double>("n90hitMin",0);
178 >   genPtMin_ = iConfig.getUntrackedParameter<double>("genPtMin",20);
179 >   emfMin_ = iConfig.getUntrackedParameter<double>("emfMin",0.01);
180 >   n90Min_ = iConfig.getUntrackedParameter<double>("n90Min",1);
181 >   n90hitMin_ = iConfig.getUntrackedParameter<double>("n90hitMin",1);
182  
183     filterJets_ = iConfig.getUntrackedParameter<bool>("filterJets",true);
184 <   diJetsOnly_ = iConfig.getUntrackedParameter<bool>("diJetsOnly",true);
185 <   matchDiJets_ = iConfig.getUntrackedParameter<bool>("usePat",true);
184 >   diJetsOnly_ = iConfig.getUntrackedParameter<bool>("diJetsOnly",false);
185 >   matchDiJets_ = iConfig.getUntrackedParameter<bool>("matchDiJets",false);
186 >   matchPatGen_ = iConfig.getUntrackedParameter<bool>("matchPatGen",false);
187 >
188 >   matchNew_ = iConfig.getUntrackedParameter<bool>("matchNew",false);
189 >
190     usePat_ = iConfig.getUntrackedParameter<bool>("usePat",true);
191     doMC_ = iConfig.getUntrackedParameter<bool>("doMC",true);
174 –   jetTag_ = iConfig.getUntrackedParameter<edm::InputTag>("src",edm::InputTag("selectedPatJets"));
192  
193 +   sortJets_ = iConfig.getUntrackedParameter<bool>("sortJets",true);
194 +   correctJets_ = iConfig.getUntrackedParameter<bool>("correctJets",false);
195 +
196 +   jetTag_ = iConfig.getUntrackedParameter<edm::InputTag>("src",edm::InputTag("selectedPatJets"));
197 +   matchTag_ = iConfig.getUntrackedParameter<edm::InputTag>("match",edm::InputTag("selectedPatJets"));
198 +   matchTags_ = iConfig.getUntrackedParameter<std::vector<edm::InputTag> >("matches",std::vector<edm::InputTag>(0));
199   }
200  
201  
# Line 197 | Line 220 | HiJetResponseAnalyzer::analyze(const edm
220  
221     iEvent.getByLabel(jetTag_,jets);
222     if(usePat_)iEvent.getByLabel(jetTag_,patjets);
223 +   std::vector<JRAV> jraV;
224  
201 –   jra_.nref = 0;
225     for(unsigned int j = 0 ; j < jets->size(); ++j){
226 <      
227 <      //      const pat::Jet& jet = (*jets)[j];
228 <      if(filterJets_ && !selectJet(j)) continue;
229 <
230 <      const reco::Jet& jet = (*jets)[j];      
231 <
232 <      jra_.jtpt[jra_.nref] = jet.pt();
233 <      jra_.jteta[jra_.nref] = jet.eta();
234 <      jra_.jtphi[jra_.nref] = jet.phi();
235 <      jra_.jtcorpt[jra_.nref] = jet.pt();
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 >     jv.index = j;
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 >       }else{
246 >         jv.refpt = -99;
247 >         jv.refeta = -99;
248 >         jv.refphi = -99;
249 >       }
250 >     }
251  
252 <      if(usePat_){
253 <         const pat::Jet& patjet = (*patjets)[j];
252 >     jraV.push_back(jv);
253 >   }
254  
255 <         jra_.jtpt[jra_.nref] = patjet.correctedJet("raw").pt();
256 <         jra_.jtcorpt[jra_.nref] = patjet.pt();
257 <        
220 <         if(doMC_ && patjet.genJet() != 0){        
221 <            if(patjet.genJet()->pt() < genPtMin_) continue;
222 <            jra_.refpt[jra_.nref] = patjet.genJet()->pt();
223 <            jra_.refeta[jra_.nref] = patjet.genJet()->eta();
224 <            jra_.refphi[jra_.nref] = patjet.genJet()->phi();
225 <                            
226 <         }else{
227 <            jra_.refpt[jra_.nref] = -99;
228 <            jra_.refeta[jra_.nref] = -99;
229 <            jra_.refphi[jra_.nref] = -99;
230 <         }
231 <      }
232 <      
233 <      jra_.nref++;
234 <      
255 >   if(sortJets_){
256 >     if(usePat_ || correctJets_) std::sort(jraV.begin(),jraV.end(),compareCorPt);
257 >     else std::sort(jraV.begin(),jraV.end(),comparePt);
258     }
259  
260 +   for(unsigned int i = 0; i < jraV.size(); ++i){
261 +     JRAV& jv = jraV[i];
262 +     const reco::Jet& jet = (*jets)[jv.index];
263 +    
264 +     if(matchNew_){
265 +       for(unsigned int im = 0; im < matchTags_.size(); ++im){
266 +         iEvent.getByLabel(matchTags_[im],matchedJets);
267 +         for(unsigned int m = 0 ; m < matchedJets->size(); ++m){
268 +           const reco::Jet& match = (*matchedJets)[m];
269 +           double dr = reco::deltaR(jet.eta(),jet.phi(),match.eta(),match.phi());
270 +           if(dr < matchR_){
271 +             jraMatch_[im].jtcorpt[i] = -99;
272 +             jraMatch_[im].jtpt[i] = match.pt();
273 +             jraMatch_[im].jteta[i] = match.eta();
274 +             jraMatch_[im].jtphi[i] = match.phi();
275 +           }
276 +         }
277 +       }
278 +     }
279 +    
280 +     jra_.jtpt[i] = jv.jtpt;
281 +     jra_.jteta[i] = jv.jteta;
282 +     jra_.jtphi[i] = jv.jtphi;
283 +     jra_.jtcorpt[i] = jv.jtcorpt;
284 +     jra_.refpt[i] = jv.refpt;
285 +     jra_.refeta[i] = jv.refeta;
286 +     jra_.refphi[i] = jv.refphi;
287 +   }
288 +   jra_.nref = jraV.size();
289 +  
290     t->Fill();
291  
292   }
293  
241 –
294   // ------------ method called once each job just before starting event loop  ------------
243 – void
244 – HiJetResponseAnalyzer::beginJob()
245 – {
295  
296 <   t= fs->make<TTree>("t","Jet Response Analyzer");
297 <   t->Branch("b",&jra_.b,"b/F");
298 <   t->Branch("hf",&jra_.hf,"hf/F");
299 <   t->Branch("nref",&jra_.nref,"nref/I");
300 <   t->Branch("jtpt",jra_.jtpt,"jtpt[nref]/F");
301 <   t->Branch("jtcorpt",jra_.jtcorpt,"jtcorpt[nref]/F");
302 <   t->Branch("refpt",jra_.refpt,"refpt[nref]/F");
303 <   t->Branch("jteta",jra_.jteta,"jteta[nref]/F");
296 >
297 > void
298 > HiJetResponseAnalyzer::beginJob(){
299 >  t= fs->make<TTree>("t","Jet Response Analyzer");
300 >  t->Branch("nref",&jra_.nref,"nref/I");
301 >  t->Branch("jtpt",jra_.jtpt,"jtpt[nref]/F");
302 >  t->Branch("jtcorpt",jra_.jtcorpt,"jtcorpt[nref]/F");
303 >  t->Branch("jteta",jra_.jteta,"jteta[nref]/F");
304 >  t->Branch("jtphi",jra_.jtphi,"jtphi[nref]/F");
305 >  t->Branch("refpt",jra_.refpt,"refpt[nref]/F");
306 >  t->Branch("refcorpt",jra_.refpt,"refcorpt[nref]/F");
307     t->Branch("refeta",jra_.refeta,"refeta[nref]/F");
256 –   t->Branch("jtphi",jra_.jtphi,"jtphi[nref]/F");
308     t->Branch("refphi",jra_.refphi,"refphi[nref]/F");
309     t->Branch("weight",&jra_.weight,"weight/F");
310 <   t->Branch("bin",&jra_.bin,"bin/I");
311 <
312 <
313 <
310 >   for(unsigned int im = 0; im < matchTags_.size(); ++im){
311 >     JRA jrm;
312 >     jraMatch_.push_back(jrm);
313 >     t->Branch(Form("jtpt%d",im),jraMatch_[im].jtpt,Form("jtpt%d[nref]/F",im));
314 >     t->Branch(Form("jtcorpt%d",im),jraMatch_[im].jtcorpt,Form("jtcorpt%d[nref]/F",im));
315 >     t->Branch(Form("jteta%d",im),jraMatch_[im].jteta,Form("jteta%d[nref]/F",im));
316 >     t->Branch(Form("jtphi%d",im),jraMatch_[im].jtphi,Form("jtphi%d[nref]/F",im));
317 >   }
318 >  
319   }
320  
321   // ------------ method called once each job just after ending the event loop  ------------

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines