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

Comparing UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc (file contents):
Revision 1.5 by peiffer, Thu Apr 5 09:48:20 2012 UTC vs.
Revision 1.19 by peiffer, Tue Jun 12 11:43:16 2012 UTC

# Line 50 | Line 50 | NtupleWriter::NtupleWriter(const edm::Pa
50    doMuons = iConfig.getParameter<bool>("doMuons");
51    doTaus = iConfig.getParameter<bool>("doTaus");
52    doJets = iConfig.getParameter<bool>("doJets");
53 +  doJECUncertainty = iConfig.getParameter<bool>("doJECUncertainty");
54 +  doGenTopJets = iConfig.getParameter<bool>("doGenTopJets");  
55    doPhotons = iConfig.getParameter<bool>("doPhotons");
56    doMET = iConfig.getParameter<bool>("doMET");
57    doGenInfo = iConfig.getParameter<bool>("doGenInfo");
58 +  doAllGenParticles = iConfig.getParameter<bool>("doAllGenParticles");
59 +  doLumiInfo = iConfig.getParameter<bool>("doLumiInfo");
60    doPV = iConfig.getParameter<bool>("doPV");
61    doTopJets = iConfig.getParameter<bool>("doTopJets");
62 +  doTrigger = iConfig.getParameter<bool>("doTrigger");
63  
64    // initialization of tree variables
65  
# Line 62 | Line 67 | NtupleWriter::NtupleWriter(const edm::Pa
67    tr->Branch("event",&event);
68    tr->Branch("luminosityBlock",&luminosityBlock);
69    tr->Branch("isRealData",&isRealData);
70 <  tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
71 <  tr->Branch("beamspot_x0",&beamspot_x0);
67 <  tr->Branch("beamspot_y0",&beamspot_y0);
68 <  tr->Branch("beamspot_z0",&beamspot_z0);
70 >  tr->Branch("rho",&rho);
71 >  rho_source = iConfig.getParameter<edm::InputTag>("rho_source");
72  
73 +  //tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
74 +  if(doLumiInfo){
75 +    tr->Branch("intgRecLumi",&intgRecLumi);
76 +    tr->Branch("intgDelLumi",&intgDelLumi);
77 +  }
78 +  if(doPV){
79 +    tr->Branch("beamspot_x0",&beamspot_x0);
80 +    tr->Branch("beamspot_y0",&beamspot_y0);
81 +    tr->Branch("beamspot_z0",&beamspot_z0);
82 +  }
83    if(doElectrons){
84      electron_sources = iConfig.getParameter<std::vector<std::string> >("electron_sources");
85      for(size_t j=0; j< electron_sources.size(); ++j){  
# Line 103 | Line 116 | NtupleWriter::NtupleWriter(const edm::Pa
116        tr->Branch( topjet_sources[j].c_str(), "std::vector<TopJet>", &topjets[j]);
117      }
118    }
119 +  if(doGenTopJets){
120 +    gentopjet_sources = iConfig.getParameter<std::vector<std::string> >("gentopjet_sources");
121 +    gentopjet_ptmin = iConfig.getParameter<double> ("gentopjet_ptmin");
122 +    gentopjet_etamax = iConfig.getParameter<double> ("gentopjet_etamax");
123 +    for(size_t j=0; j< gentopjet_sources.size(); ++j){  
124 +      tr->Branch( gentopjet_sources[j].c_str(), "std::vector<TopJet>", &gentopjets[j]);
125 +    }
126 +  }
127    if(doPhotons){
128      photon_sources = iConfig.getParameter<std::vector<std::string> >("photon_sources");
129      for(size_t j=0; j< photon_sources.size(); ++j){  
# Line 125 | Line 146 | NtupleWriter::NtupleWriter(const edm::Pa
146      tr->Branch("genInfo","GenInfo",&genInfo);
147      tr->Branch("GenParticles","std::vector<GenParticle>", &genps);
148    }
149 <
150 <  trigger_prefixes = iConfig.getParameter<std::vector<std::string> >("trigger_prefixes");
151 <  //tr->Branch("triggerResults","std::map<std::string, bool>",&triggerResults);
152 <  tr->Branch("triggerNames", "std::vector<std::string>", &triggerNames);  
153 <  tr->Branch("triggerResults", "std::vector<bool>", &triggerResults);
154 <  tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
155 <  tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
156 <
149 >  if(doTrigger){
150 >    trigger_prefixes = iConfig.getParameter<std::vector<std::string> >("trigger_prefixes");
151 >    //tr->Branch("triggerResults","std::map<std::string, bool>",&triggerResults);
152 >    tr->Branch("triggerNames", "std::vector<std::string>", &triggerNames);  
153 >    tr->Branch("triggerResults", "std::vector<bool>", &triggerResults);
154 >    //tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
155 >    //tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
156 >  }
157    newrun = true;
158   }
159  
# Line 164 | Line 185 | NtupleWriter::analyze(const edm::Event&
185     luminosityBlock = iEvent.luminosityBlock();
186     isRealData      = iEvent.isRealData();
187  
188 <   if(isRealData){
189 <     edm::Handle<bool> bool_handle;
190 <     iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
191 <     HBHENoiseFilterResult = *(bool_handle.product());
192 <   }
193 <   else HBHENoiseFilterResult = false;
188 >   edm::Handle<double> m_rho;
189 >   iEvent.getByLabel(rho_source,m_rho);
190 >   rho=*m_rho;
191 >
192 > //    if(isRealData){
193 > //      edm::Handle<bool> bool_handle;
194 > //      iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
195 > //      HBHENoiseFilterResult = *(bool_handle.product());
196 > //    }
197 > //    else HBHENoiseFilterResult = false;
198  
199     // ------------- primary vertices and beamspot  -------------
200  
# Line 185 | Line 210 | NtupleWriter::analyze(const edm::Event&
210           reco::Vertex reco_pv = reco_pvs[i];
211  
212           PrimaryVertex pv;
213 <         pv.x =  reco_pv.x();
214 <         pv.y =  reco_pv.y();
215 <         pv.z =  reco_pv.z();
216 <         pv.nTracks =  reco_pv.nTracks();
217 <         //pv.isValid =  reco_pv.isValid();
218 <         pv.chi2 =  reco_pv.chi2();
219 <         pv.ndof =  reco_pv.ndof();      
213 >         pv.set_x( reco_pv.x());
214 >         pv.set_y( reco_pv.y());
215 >         pv.set_z( reco_pv.z());
216 >         pv.set_nTracks( reco_pv.nTracks());
217 >         //pv.set_isValid( reco_pv.isValid());
218 >         pv.set_chi2( reco_pv.chi2());
219 >         pv.set_ndof( reco_pv.ndof());  
220  
221           pvs[j].push_back(pv);
222         }
223       }
224 +      
225 +     edm::Handle<reco::BeamSpot> beamSpot;
226 +     iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
227 +     const reco::BeamSpot & bsp = *beamSpot;
228 +    
229 +     beamspot_x0 = bsp.x0();
230 +     beamspot_y0 = bsp.y0();
231 +     beamspot_z0 = bsp.z0();
232     }
233 +
234 + // ------------- generator info -------------
235    
236 <   edm::Handle<reco::BeamSpot> beamSpot;
237 <   iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
238 <   const reco::BeamSpot & bsp = *beamSpot;
239 <  
240 <   beamspot_x0 = bsp.x0();
241 <   beamspot_y0 = bsp.y0();
242 <   beamspot_z0 = bsp.z0();
236 >   if(doGenInfo){
237 >     genInfo.clear_weights();
238 >     genInfo.clear_binningValues();
239 >     genps.clear();
240 >
241 >     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
242 >     iEvent.getByLabel("generator", genEventInfoProduct);
243 >     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
244 >  
245 >     for(unsigned int k=0; k<genEventInfo.binningValues().size();++k){
246 >       genInfo.add_binningValue(genEventInfo.binningValues().at(k));
247 >     }
248 >     for(unsigned int k=0; k<genEventInfo.weights().size();++k){
249 >       genInfo.add_weight(genEventInfo.weights().at(k));
250 >     }
251 >     genInfo.set_alphaQCD(genEventInfo.alphaQCD());
252 >     genInfo.set_alphaQED(genEventInfo.alphaQED());
253 >     genInfo.set_qScale(genEventInfo.qScale());
254 >    
255 >     const gen::PdfInfo* pdf = genEventInfo.pdf();
256 >     if(pdf){
257 >       genInfo.set_pdf_id1(pdf->id.first);
258 >       genInfo.set_pdf_id2(pdf->id.second);
259 >       genInfo.set_pdf_x1(pdf->x.first);
260 >       genInfo.set_pdf_x2(pdf->x.second);
261 >       genInfo.set_pdf_scalePDF(pdf->scalePDF);
262 >       genInfo.set_pdf_xPDF1(pdf->xPDF.first);
263 >       genInfo.set_pdf_xPDF2(pdf->xPDF.second);
264 >     }
265 >     else{
266 >       genInfo.set_pdf_id1(-999);
267 >       genInfo.set_pdf_id2(-999);
268 >       genInfo.set_pdf_x1(-999);
269 >       genInfo.set_pdf_x2(-999);
270 >       genInfo.set_pdf_scalePDF(-999);
271 >       genInfo.set_pdf_xPDF1(-999);
272 >       genInfo.set_pdf_xPDF2(-999);
273 >     }
274 >
275 >     edm::Handle<std::vector<PileupSummaryInfo> > pus;
276 >     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
277 >     genInfo.set_pileup_NumInteractions_intime(0);
278 >     genInfo.set_pileup_NumInteractions_ootbefore(0);
279 >     genInfo.set_pileup_NumInteractions_ootafter(0);
280 >     if(pus.isValid()){
281 >       genInfo.set_pileup_TrueNumInteractions ( (float) pus->at(0).getTrueNumInteractions());
282 >       for(size_t i=0; i<pus->size(); ++i){
283 >         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
284 >           genInfo.set_pileup_NumInteractions_intime( genInfo.pileup_NumInteractions_intime() + pus->at(i).getPU_NumInteractions());
285 >         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
286 >           genInfo.set_pileup_NumInteractions_ootbefore( genInfo.pileup_NumInteractions_ootbefore() + pus->at(i).getPU_NumInteractions());
287 >         }
288 >         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
289 >           genInfo.set_pileup_NumInteractions_ootafter( genInfo.pileup_NumInteractions_ootafter() + pus->at(i).getPU_NumInteractions());
290 >         }
291 >       }
292 >     }
293 >
294 >     edm::Handle<reco::GenParticleCollection> genPartColl;
295 >     iEvent.getByLabel(edm::InputTag("genParticles"), genPartColl);
296 >     int index=-1;
297 >     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
298 >       index++;
299 >      
300 >       //write out only top quarks and status 3 particles (works fine only for MadGraph)
301 >       if(abs(iter->pdgId())==6 || iter->status()==3 || doAllGenParticles){
302 >         GenParticle genp;
303 >         genp.set_charge(iter->charge());
304 >         genp.set_pt(iter->p4().pt());
305 >         genp.set_eta(iter->p4().eta());
306 >         genp.set_phi(iter->p4().phi());
307 >         genp.set_energy(iter->p4().E());
308 >         genp.set_index(index);
309 >         genp.set_status( iter->status());
310 >         genp.set_pdgId( iter->pdgId());
311 >
312 >         genp.set_mother1(-1);
313 >         genp.set_mother2(-1);
314 >         genp.set_daughter1(-1);
315 >         genp.set_daughter2(-1);
316 >        
317 >         int nm=iter->numberOfMothers();
318 >         int nd=iter->numberOfDaughters();
319 >
320 >        
321 >         if (nm>0) genp.set_mother1( iter->motherRef(0).key());
322 >         if (nm>1) genp.set_mother2( iter->motherRef(1).key());
323 >         if (nd>0) genp.set_daughter1( iter->daughterRef(0).key());
324 >         if (nd>1) genp.set_daughter2( iter->daughterRef(1).key());
325 >
326 >         genps.push_back(genp);
327 >       }
328 >     }
329 >   }
330  
331     // ------------- electrons -------------  
332     if(doElectrons){
333 +
334 +     edm::Handle<reco::ConversionCollection> hConversions;
335 +     iEvent.getByLabel("allConversions", hConversions);
336 +
337 +     edm::Handle<reco::BeamSpot> beamSpot;
338 +     iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
339 +     const reco::BeamSpot & bsp = *beamSpot;
340 +
341       for(size_t j=0; j< electron_sources.size(); ++j){
342         eles[j].clear();
343         edm::Handle< std::vector<pat::Electron> > ele_handle;
# Line 218 | Line 348 | NtupleWriter::analyze(const edm::Event&
348           pat::Electron pat_ele = pat_electrons[i];
349           Electron ele;
350          
351 <         ele.charge =  pat_ele.charge();
352 <         ele.pt =  pat_ele.pt();
353 <         ele.eta =  pat_ele.eta();
354 <         ele.phi =  pat_ele.phi();
355 <         ele.energy =  pat_ele.energy();
356 <         ele.vertex_x = pat_ele.vertex().x();
357 <         ele.vertex_y = pat_ele.vertex().y();
358 <         ele.vertex_z = pat_ele.vertex().z();
359 <         ele.supercluster_eta = pat_ele.superCluster()->eta();
360 <         ele.supercluster_phi = pat_ele.superCluster()->phi();
361 <         ele.dB = pat_ele.dB();
362 <         //ele.particleIso = pat_ele.particleIso();
363 <         ele.neutralHadronIso = pat_ele.neutralHadronIso();
364 <         ele.chargedHadronIso = pat_ele.chargedHadronIso();
365 <         ele.trackIso = pat_ele.trackIso();
366 <         ele.puChargedHadronIso = pat_ele.puChargedHadronIso();
367 <         ele.gsfTrack_trackerExpectedHitsInner_numberOfLostHits = pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits();
368 <         ele.gsfTrack_px= pat_ele.gsfTrack()->px();
369 <         ele.gsfTrack_py= pat_ele.gsfTrack()->py();
370 <         ele.gsfTrack_pz= pat_ele.gsfTrack()->pz();
371 <         ele.gsfTrack_vx= pat_ele.gsfTrack()->vx();
372 <         ele.gsfTrack_vy= pat_ele.gsfTrack()->vy();
373 <         ele.gsfTrack_vz= pat_ele.gsfTrack()->vz();
351 >         ele.set_charge( pat_ele.charge());
352 >         ele.set_pt( pat_ele.pt());
353 >         ele.set_eta( pat_ele.eta());
354 >         ele.set_phi( pat_ele.phi());
355 >         ele.set_energy( pat_ele.energy());
356 >         ele.set_vertex_x(pat_ele.vertex().x());
357 >         ele.set_vertex_y(pat_ele.vertex().y());
358 >         ele.set_vertex_z(pat_ele.vertex().z());
359 >         ele.set_supercluster_eta(pat_ele.superCluster()->eta());
360 >         ele.set_supercluster_phi(pat_ele.superCluster()->phi());
361 >         ele.set_dB(pat_ele.dB());
362 >         //ele.set_particleIso(pat_ele.particleIso());
363 >         ele.set_neutralHadronIso(pat_ele.neutralHadronIso());
364 >         ele.set_chargedHadronIso(pat_ele.chargedHadronIso());
365 >         ele.set_trackIso(pat_ele.trackIso());
366 >         ele.set_photonIso(pat_ele.photonIso());
367 >         ele.set_puChargedHadronIso(pat_ele.puChargedHadronIso());
368 >         ele.set_gsfTrack_trackerExpectedHitsInner_numberOfLostHits(pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits());
369 >         ele.set_gsfTrack_px( pat_ele.gsfTrack()->px());
370 >         ele.set_gsfTrack_py( pat_ele.gsfTrack()->py());
371 >         ele.set_gsfTrack_pz( pat_ele.gsfTrack()->pz());
372 >         ele.set_gsfTrack_vx( pat_ele.gsfTrack()->vx());
373 >         ele.set_gsfTrack_vy( pat_ele.gsfTrack()->vy());
374 >         ele.set_gsfTrack_vz( pat_ele.gsfTrack()->vz());
375 >         ele.set_passconversionveto(!ConversionTools::hasMatchedConversion(pat_ele,hConversions,bsp.position()));
376 >         ele.set_dEtaIn(pat_ele.deltaEtaSuperClusterTrackAtVtx());
377 >         ele.set_dPhiIn(pat_ele.deltaPhiSuperClusterTrackAtVtx());
378 >         ele.set_sigmaIEtaIEta(pat_ele.sigmaIetaIeta());
379 >         ele.set_HoverE(pat_ele.hadronicOverEm());
380 >         ele.set_fbrem(pat_ele.fbrem());
381 >         ele.set_EoverPIn(pat_ele.eSuperClusterOverP());
382 >         ele.set_EcalEnergy(pat_ele.ecalEnergy());
383 >         ele.set_mvaTrigV0(pat_ele.electronID("mvaTrigV0"));
384 >         ele.set_mvaNonTrigV0(pat_ele.electronID("mvaNonTrigV0"));
385 >
386           eles[j].push_back(ele);
387         }
388       }
# Line 259 | Line 401 | NtupleWriter::analyze(const edm::Event&
401           pat::Muon pat_mu = pat_muons[i];
402  
403           Muon mu;
404 <         mu.charge =  pat_mu.charge();
405 <         mu.pt =  pat_mu.pt();
406 <         mu.eta =  pat_mu.eta();
407 <         mu.phi =  pat_mu.phi();
408 <         mu.energy =  pat_mu.energy();
409 <         mu.vertex_x = pat_mu.vertex().x();
410 <         mu.vertex_y = pat_mu.vertex().y();
411 <         mu.vertex_z = pat_mu.vertex().z();
412 <         mu.dB = pat_mu.dB();
413 <         //mu.particleIso = pat_mu.particleIso();
414 <         mu.neutralHadronIso = pat_mu.neutralHadronIso();
415 <         mu.chargedHadronIso = pat_mu.chargedHadronIso();
416 <         mu.trackIso = pat_mu.trackIso();
417 <         mu.puChargedHadronIso = pat_mu.puChargedHadronIso();
418 <         mu.isGlobalMuon = pat_mu.isGlobalMuon();
419 <         mu.isStandAloneMuon = pat_mu.isStandAloneMuon();
420 <         mu.isTrackerMuon = pat_mu.isTrackerMuon();
421 <         mu.numberOfMatchedStations = pat_mu.numberOfMatchedStations();
404 >         mu.set_charge( pat_mu.charge());
405 >         mu.set_pt( pat_mu.pt());
406 >         mu.set_eta( pat_mu.eta());
407 >         mu.set_phi( pat_mu.phi());
408 >         mu.set_energy( pat_mu.energy());
409 >         mu.set_vertex_x ( pat_mu.vertex().x());
410 >         mu.set_vertex_y ( pat_mu.vertex().y());
411 >         mu.set_vertex_z ( pat_mu.vertex().z());
412 >         mu.set_dB ( pat_mu.dB());
413 >         //mu.particleIso ( pat_mu.particleIso());
414 >         mu.set_neutralHadronIso ( pat_mu.neutralHadronIso());
415 >         mu.set_chargedHadronIso ( pat_mu.chargedHadronIso());
416 >         mu.set_trackIso ( pat_mu.trackIso());
417 >         mu.set_photonIso ( pat_mu.photonIso());
418 >         mu.set_puChargedHadronIso ( pat_mu.puChargedHadronIso());
419 >         mu.set_isGlobalMuon ( pat_mu.isGlobalMuon());
420 >         mu.set_isStandAloneMuon ( pat_mu.isStandAloneMuon());
421 >         mu.set_isTrackerMuon ( pat_mu.isTrackerMuon());
422 >         mu.set_numberOfMatchedStations ( pat_mu.numberOfMatchedStations());
423           reco::TrackRef globalTrack = pat_mu.globalTrack();
424           if(!globalTrack.isNull()){
425 <           mu.globalTrack_chi2 = globalTrack->chi2();
426 <           mu.globalTrack_ndof = globalTrack->ndof();
427 <           mu.globalTrack_d0 = globalTrack->d0();        
428 <           mu.globalTrack_d0Error = globalTrack->d0Error();
429 <           mu.globalTrack_numberOfValidHits = globalTrack->numberOfValidHits();
430 <           mu.globalTrack_numberOfLostHits = globalTrack->numberOfLostHits();
425 >           mu.set_globalTrack_chi2 ( globalTrack->chi2());
426 >           mu.set_globalTrack_ndof ( globalTrack->ndof());
427 >           mu.set_globalTrack_d0 ( globalTrack->d0());  
428 >           mu.set_globalTrack_d0Error ( globalTrack->d0Error());
429 >           mu.set_globalTrack_numberOfValidHits ( globalTrack->numberOfValidHits());
430 >           mu.set_globalTrack_numberOfLostHits ( globalTrack->numberOfLostHits());
431           }
432           else{
433 <           mu.globalTrack_chi2 = 0;
434 <           mu.globalTrack_ndof = 0;
435 <           mu.globalTrack_d0 = 0;
436 <           mu.globalTrack_d0Error = 0;
437 <           mu.globalTrack_numberOfValidHits = 0;
438 <           mu.globalTrack_numberOfLostHits = 0;
433 >           mu.set_globalTrack_chi2 ( 0);
434 >           mu.set_globalTrack_ndof ( 0);
435 >           mu.set_globalTrack_d0 ( 0);
436 >           mu.set_globalTrack_d0Error ( 0);
437 >           mu.set_globalTrack_numberOfValidHits ( 0);
438 >           mu.set_globalTrack_numberOfLostHits ( 0);
439           }
440           reco::TrackRef innerTrack = pat_mu.innerTrack();
441           if(!innerTrack.isNull()){
442 <           mu.innerTrack_chi2 = innerTrack->chi2();
443 <           mu.innerTrack_ndof = innerTrack->ndof();
444 <           mu.innerTrack_d0 = innerTrack->d0();  
445 <           mu.innerTrack_d0Error = innerTrack->d0Error();
446 <           mu.innerTrack_numberOfValidHits = innerTrack->numberOfValidHits();
447 <           mu.innerTrack_numberOfLostHits = innerTrack->numberOfLostHits();
442 >           mu.set_innerTrack_chi2 ( innerTrack->chi2());
443 >           mu.set_innerTrack_ndof ( innerTrack->ndof());
444 >           mu.set_innerTrack_d0 ( innerTrack->d0());    
445 >           mu.set_innerTrack_d0Error ( innerTrack->d0Error());
446 >           mu.set_innerTrack_numberOfValidHits ( innerTrack->numberOfValidHits());
447 >           mu.set_innerTrack_numberOfLostHits ( innerTrack->numberOfLostHits());
448 >           mu.set_innerTrack_trackerLayersWithMeasurement ( innerTrack->hitPattern().trackerLayersWithMeasurement());
449 >           mu.set_innerTrack_numberOfValidPixelHits ( innerTrack->hitPattern().numberOfValidPixelHits());
450           }
451           else{
452 <           mu.innerTrack_chi2 = 0;
453 <           mu.innerTrack_ndof = 0;
454 <           mu.innerTrack_d0 = 0;
455 <           mu.innerTrack_d0Error = 0;
456 <           mu.innerTrack_numberOfValidHits = 0;
457 <           mu.innerTrack_numberOfLostHits = 0;
452 >           mu.set_innerTrack_chi2 ( 0);
453 >           mu.set_innerTrack_ndof ( 0);
454 >           mu.set_innerTrack_d0 ( 0);
455 >           mu.set_innerTrack_d0Error ( 0);
456 >           mu.set_innerTrack_numberOfValidHits ( 0);
457 >           mu.set_innerTrack_numberOfLostHits ( 0);
458 >           mu.set_innerTrack_trackerLayersWithMeasurement ( 0);
459 >           mu.set_innerTrack_numberOfValidPixelHits ( 0);
460           }
461           reco::TrackRef outerTrack = pat_mu.outerTrack();
462           if(!outerTrack.isNull()){
463 <           mu.outerTrack_chi2 = outerTrack->chi2();
464 <           mu.outerTrack_ndof = outerTrack->ndof();
465 <           mu.outerTrack_d0 = outerTrack->d0();  
466 <           mu.outerTrack_d0Error = outerTrack->d0Error();
467 <           mu.outerTrack_numberOfValidHits = outerTrack->numberOfValidHits();
468 <           mu.outerTrack_numberOfLostHits = outerTrack->numberOfLostHits();
463 >           mu.set_outerTrack_chi2 ( outerTrack->chi2());
464 >           mu.set_outerTrack_ndof ( outerTrack->ndof());
465 >           mu.set_outerTrack_d0 ( outerTrack->d0());    
466 >           mu.set_outerTrack_d0Error ( outerTrack->d0Error());
467 >           mu.set_outerTrack_numberOfValidHits ( outerTrack->numberOfValidHits());
468 >           mu.set_outerTrack_numberOfLostHits ( outerTrack->numberOfLostHits());
469           }
470           else{
471 <           mu.outerTrack_chi2 = 0;
472 <           mu.outerTrack_ndof = 0;
473 <           mu.outerTrack_d0 = 0;
474 <           mu.outerTrack_d0Error = 0;
475 <           mu.outerTrack_numberOfValidHits = 0;
476 <           mu.outerTrack_numberOfLostHits = 0;
471 >           mu.set_outerTrack_chi2 ( 0);
472 >           mu.set_outerTrack_ndof ( 0);
473 >           mu.set_outerTrack_d0 ( 0);
474 >           mu.set_outerTrack_d0Error ( 0);
475 >           mu.set_outerTrack_numberOfValidHits ( 0);
476 >           mu.set_outerTrack_numberOfLostHits ( 0);
477           }
478  
479           mus[j].push_back(mu);
# Line 350 | Line 497 | NtupleWriter::analyze(const edm::Event&
497           if(fabs(pat_tau.eta()) > tau_etamax) continue;
498  
499           Tau tau;
500 <         tau.charge =  pat_tau.charge();
501 <         tau.pt =  pat_tau.pt();
502 <         tau.eta =  pat_tau.eta();
503 <         tau.phi =  pat_tau.phi();
504 <         tau.energy =  pat_tau.energy();
505 <
500 >         tau.set_charge( pat_tau.charge());
501 >         tau.set_pt( pat_tau.pt());
502 >         tau.set_eta( pat_tau.eta());
503 >         tau.set_phi( pat_tau.phi());
504 >         tau.set_energy( pat_tau.energy());
505 >         tau.set_decayModeFinding ( pat_tau.tauID("decayModeFinding")>0.5);
506 >         tau.set_byVLooseCombinedIsolationDeltaBetaCorr  ( pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5);
507 >         tau.set_byLooseCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5);
508 >         tau.set_byMediumCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5);
509 >         tau.set_byTightCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5);
510 >         tau.set_againstElectronLoose  ( pat_tau.tauID("againstElectronLoose")>0.5);
511 >         tau.set_againstElectronMedium ( pat_tau.tauID("againstElectronMedium")>0.5);
512 >         tau.set_againstElectronTight ( pat_tau.tauID("againstElectronTight")>0.5);
513 >         tau.set_againstElectronMVA  ( pat_tau.tauID("againstElectronMVA")>0.5);
514 >         tau.set_againstMuonLoose ( pat_tau.tauID("againstMuonLoose")>0.5);
515 >         tau.set_againstMuonMedium ( pat_tau.tauID("againstMuonMedium")>0.5);
516 >         tau.set_againstMuonTight ( pat_tau.tauID("againstMuonTight")>0.5);
517 >
518 > //       reco::PFCandidateRef leadPFCand = pat_tau.leadPFCand();
519 > //       if(!leadPFCand.isNull()){
520 > //         tau.set_leadPFCand_px ( leadPFCand->px());
521 > //         tau.set_leadPFCand_py ( leadPFCand->py());
522 > //         tau.set_leadPFCand_pz ( leadPFCand->pz());
523 > //       }
524 > //       else{
525 > //         tau.set_leadPFCand_px ( 0);
526 > //         tau.set_leadPFCand_py ( 0);
527 > //         tau.set_leadPFCand_pz ( 0);
528 > //       }
529           taus[j].push_back(tau);
530         }
531       }
# Line 382 | Line 552 | NtupleWriter::analyze(const edm::Event&
552   //       }
553  
554           Jet jet;
555 <         jet.charge = pat_jet.charge();
556 <         jet.pt = pat_jet.pt();
557 <         jet.eta = pat_jet.eta();
558 <         jet.phi = pat_jet.phi();
559 <         jet.energy = pat_jet.energy();
560 <         jet.numberOfDaughters =pat_jet.numberOfDaughters();
555 >         jet.set_charge(pat_jet.charge());
556 >         jet.set_pt(pat_jet.pt());
557 >         jet.set_eta(pat_jet.eta());
558 >         jet.set_phi(pat_jet.phi());
559 >         jet.set_energy(pat_jet.energy());
560 >         jet.set_numberOfDaughters (pat_jet.numberOfDaughters());
561           const reco::TrackRefVector&  jettracks = pat_jet.associatedTracks();
562 <         jet.nTracks = jettracks.size();
563 <         jet.jetArea = pat_jet.jetArea();
564 <         jet.pileup = pat_jet.pileup();
562 >         jet.set_nTracks ( jettracks.size());
563 >         jet.set_jetArea(pat_jet.jetArea());
564 >         jet.set_pileup(pat_jet.pileup());
565           if(pat_jet.isPFJet()){
566 <           jet.neutralEmEnergyFraction =pat_jet.neutralEmEnergyFraction();
567 <           jet.neutralHadronEnergyFraction =pat_jet.neutralHadronEnergyFraction();
568 <           jet.chargedEmEnergyFraction =pat_jet.chargedEmEnergyFraction();
569 <           jet.chargedHadronEnergyFraction =pat_jet.chargedHadronEnergyFraction();
570 <           jet.muonEnergyFraction =pat_jet.muonEnergyFraction();
571 <           jet.photonEnergyFraction =pat_jet.photonEnergyFraction();
572 <           jet.chargedMultiplicity =pat_jet.chargedMultiplicity();
573 <           jet.neutralMultiplicity =pat_jet.neutralMultiplicity();
574 <           jet.muonMultiplicity =pat_jet.muonMultiplicity();
575 <           jet.electronMultiplicity =pat_jet.electronMultiplicity();
576 <           jet.photonMultiplicity =pat_jet.photonMultiplicity();
577 <         }
578 <
579 <         jecUnc->setJetEta(pat_jet.eta());
580 <         jecUnc->setJetPt(pat_jet.pt());
581 <         jet.JEC_uncertainty = jecUnc->getUncertainty(true);
582 <
583 <         jet.btag_simpleSecondaryVertexHighEff=pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
584 <         jet.btag_simpleSecondaryVertexHighPur=pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
585 <         jet.btag_combinedSecondaryVertex=pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags");
586 <         jet.btag_combinedSecondaryVertexMVA=pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
587 <         jet.btag_jetBProbability=pat_jet.bDiscriminator("jetBProbabilityBJetTags");
588 <         jet.btag_jetProbability=pat_jet.bDiscriminator("jetProbabilityBJetTags");
589 <
566 >           jet.set_neutralEmEnergyFraction (pat_jet.neutralEmEnergyFraction());
567 >           jet.set_neutralHadronEnergyFraction (pat_jet.neutralHadronEnergyFraction());
568 >           jet.set_chargedEmEnergyFraction (pat_jet.chargedEmEnergyFraction());
569 >           jet.set_chargedHadronEnergyFraction (pat_jet.chargedHadronEnergyFraction());
570 >           jet.set_muonEnergyFraction (pat_jet.muonEnergyFraction());
571 >           jet.set_photonEnergyFraction (pat_jet.photonEnergyFraction());
572 >           jet.set_chargedMultiplicity (pat_jet.chargedMultiplicity());
573 >           jet.set_neutralMultiplicity (pat_jet.neutralMultiplicity());
574 >           jet.set_muonMultiplicity (pat_jet.muonMultiplicity());
575 >           jet.set_electronMultiplicity (pat_jet.electronMultiplicity());
576 >           jet.set_photonMultiplicity (pat_jet.photonMultiplicity());
577 >         }
578 >         if(doJECUncertainty){
579 >           jecUnc->setJetEta(pat_jet.eta());
580 >           jecUnc->setJetPt(pat_jet.pt());
581 >           jet.set_JEC_uncertainty(jecUnc->getUncertainty(true));
582 >         }
583 >         jet.set_JEC_factor_raw(pat_jet.jecFactor("Uncorrected"));
584 >        
585 >         jet.set_btag_simpleSecondaryVertexHighEff(pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
586 >         jet.set_btag_simpleSecondaryVertexHighPur(pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
587 >         jet.set_btag_combinedSecondaryVertex(pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags"));
588 >         jet.set_btag_combinedSecondaryVertexMVA(pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
589 >         jet.set_btag_jetBProbability(pat_jet.bDiscriminator("jetBProbabilityBJetTags"));
590 >         jet.set_btag_jetProbability(pat_jet.bDiscriminator("jetProbabilityBJetTags"));
591 >
592 >         const reco::GenJet *genj = pat_jet.genJet();
593 >         if(genj){
594 >           jet.set_genjet_pt(genj->pt());
595 >           jet.set_genjet_eta(genj->eta());
596 >           jet.set_genjet_phi(genj->phi());
597 >           jet.set_genjet_energy(genj->energy());
598 >           if(doAllGenParticles){
599 >             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
600 >             for(unsigned int l = 0; l<jetgenps.size(); ++l){
601 >               for(unsigned int k=0; k< genps.size(); ++k){
602 >                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
603 >                   jet.add_genparticles_index(genps[k].index());
604 >                   break;
605 >                 }
606 >               }
607 >             }
608 >             if(jet.genparticles_indices().size()!= jetgenps.size())
609 >               std::cout << "WARNING: Found only " << jet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this jet"<<std::endl;
610 >           }
611 >         }
612           jets[j].push_back(jet);
613         }
614       }
# Line 433 | Line 625 | NtupleWriter::analyze(const edm::Event&
625         iEvent.getByLabel(topjet_sources[j], pat_topjets);
626  
627         for (unsigned int i = 0; i < pat_topjets->size(); i++) {
628 +
629           const pat::Jet  pat_topjet =  * dynamic_cast<pat::Jet const *>(&pat_topjets->at(i));
630           if(pat_topjet.pt() < topjet_ptmin) continue;
631           if(fabs(pat_topjet.eta()) > topjet_etamax) continue;
632  
633           TopJet topjet;
634 <         topjet.charge = pat_topjet.charge();
635 <         topjet.pt = pat_topjet.pt();
636 <         topjet.eta = pat_topjet.eta();
637 <         topjet.phi = pat_topjet.phi();
638 <         topjet.energy = pat_topjet.energy();
639 <         topjet.numberOfDaughters =pat_topjet.numberOfDaughters();
634 >         topjet.set_charge(pat_topjet.charge());
635 >         topjet.set_pt(pat_topjet.pt());
636 >         topjet.set_eta(pat_topjet.eta());
637 >         topjet.set_phi(pat_topjet.phi());
638 >         topjet.set_energy(pat_topjet.energy());
639 >         topjet.set_numberOfDaughters(pat_topjet.numberOfDaughters());
640           const reco::TrackRefVector&  topjettracks = pat_topjet.associatedTracks();
641 <         topjet.nTracks = topjettracks.size();
642 <         topjet.jetArea = pat_topjet.jetArea();
643 <         topjet.pileup = pat_topjet.pileup();
644 < //       topjet.neutralEmEnergyFraction =pat_topjet.neutralEmEnergyFraction();
645 < //       topjet.neutralHadronEnergyFraction =pat_topjet.neutralHadronEnergyFraction();
646 < //       topjet.chargedEmEnergyFraction =pat_topjet.chargedEmEnergyFraction();
647 < //       topjet.chargedHadronEnergyFraction =pat_topjet.chargedHadronEnergyFraction();
648 < //       topjet.muonEnergyFraction =pat_topjet.muonEnergyFraction();
649 < //       topjet.photonEnergyFraction =pat_topjet.photonEnergyFraction();
650 < //       topjet.chargedMultiplicity =pat_topjet.chargedMultiplicity();
651 < //       topjet.neutralMultiplicity =pat_topjet.neutralMultiplicity();
652 < //       topjet.muonMultiplicity =pat_topjet.muonMultiplicity();
653 < //       topjet.electronMultiplicity =pat_topjet.electronMultiplicity();
654 < //       topjet.photonMultiplicity =pat_topjet.photonMultiplicity();
655 <
656 <         jecUnc->setJetEta(pat_topjet.eta());
657 <         jecUnc->setJetPt(pat_topjet.pt());
658 <         topjet.JEC_uncertainty = jecUnc->getUncertainty(true);
659 <
660 <         topjet.btag_simpleSecondaryVertexHighEff=pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
661 <         topjet.btag_simpleSecondaryVertexHighPur=pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
662 <         topjet.btag_combinedSecondaryVertex=pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags");
663 <         topjet.btag_combinedSecondaryVertexMVA=pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
664 <         topjet.btag_jetBProbability=pat_topjet.bDiscriminator("jetBProbabilityBJetTags");
665 <         topjet.btag_jetProbability=pat_topjet.bDiscriminator("jetProbabilityBJetTags");
641 >         topjet.set_nTracks( topjettracks.size());
642 >         topjet.set_jetArea( pat_topjet.jetArea());
643 >         topjet.set_pileup( pat_topjet.pileup());
644 > //       topjet.set_neutralEmEnergyFraction(pat_topjet.neutralEmEnergyFraction());
645 > //       topjet.set_neutralHadronEnergyFraction(pat_topjet.neutralHadronEnergyFraction());
646 > //       topjet.set_chargedEmEnergyFraction(pat_topjet.chargedEmEnergyFraction());
647 > //       topjet.set_chargedHadronEnergyFraction(pat_topjet.chargedHadronEnergyFraction());
648 > //       topjet.set_muonEnergyFraction(pat_topjet.muonEnergyFraction());
649 > //       topjet.set_photonEnergyFraction(pat_topjet.photonEnergyFraction());
650 > //       topjet.set_chargedMultiplicity(pat_topjet.chargedMultiplicity());
651 > //       topjet.set_neutralMultiplicity(pat_topjet.neutralMultiplicity());
652 > //       topjet.set_muonMultiplicity(pat_topjet.muonMultiplicity());
653 > //       topjet.set_electronMultiplicity(pat_topjet.electronMultiplicity());
654 > //       topjet.set_photonMultiplicity(pat_topjet.photonMultiplicity());
655 >         if(doJECUncertainty){
656 >           jecUnc->setJetEta(pat_topjet.eta());
657 >           jecUnc->setJetPt(pat_topjet.pt());
658 >           topjet.set_JEC_uncertainty( jecUnc->getUncertainty(true));
659 >         }
660 >         topjet.set_JEC_factor_raw( pat_topjet.jecFactor("Uncorrected"));
661 >
662 >         topjet.set_btag_simpleSecondaryVertexHighEff(pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
663 >         topjet.set_btag_simpleSecondaryVertexHighPur(pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
664 >         topjet.set_btag_combinedSecondaryVertex(pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags"));
665 >         topjet.set_btag_combinedSecondaryVertexMVA(pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
666 >         topjet.set_btag_jetBProbability(pat_topjet.bDiscriminator("jetBProbabilityBJetTags"));
667 >         topjet.set_btag_jetProbability(pat_topjet.bDiscriminator("jetProbabilityBJetTags"));
668 >
669 >         /*
670 >         const reco::GenJet *genj = pat_topjet.genJet();
671 >         if(genj){
672 >           topjet.set_genjet_pt ( genj->pt());
673 >           topjet.set_genjet_eta ( genj->eta());
674 >           topjet.set_genjet_phi ( genj->phi());
675 >           topjet.set_genjet_energy ( genj->energy());
676 >           if(doAllGenParticles){
677 >             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
678 >             for(unsigned int l = 0; l<jetgenps.size(); ++l){
679 >               for(unsigned int k=0; k< genps.size(); ++k){
680 >                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
681 >                   topjet.add_genparticles_index(genps[k].index());
682 >                   break;
683 >                 }
684 >               }
685 >             }
686 >             if(topjet.genparticles_indices().size()!= jetgenps.size())
687 >               std::cout << "WARNING: Found only " << topjet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this topjet"<<std::endl;
688 >           }
689 >         }
690 >         */
691  
692           for (unsigned int k = 0; k < pat_topjet.numberOfDaughters(); k++) {
693             Particle subjet_v4;
694 <           subjet_v4.pt = pat_topjet.daughter(k)->p4().pt();
695 <           subjet_v4.eta = pat_topjet.daughter(k)->p4().eta();
696 <           subjet_v4.phi = pat_topjet.daughter(k)->p4().phi();
697 <           subjet_v4.energy = pat_topjet.daughter(k)->p4().E();
698 <           topjet.subjets.push_back(subjet_v4);
694 >           subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt());
695 >           subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta());
696 >           subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi());
697 >           subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E());
698 >           topjet.add_subjet(subjet_v4);
699           }
700           topjets[j].push_back(topjet);
701         }
702       }
703     }
704  
705 +
706 +   // ------------- generator top jets -------------
707 +   if(doGenTopJets){
708 +     for(size_t j=0; j< gentopjet_sources.size(); ++j){
709 +      
710 +       gentopjets[j].clear();
711 +
712 +       edm::Handle<reco::BasicJetCollection> reco_gentopjets;
713 +       //edm::Handle<std::vector<reco::Jet> > reco_gentopjets;
714 +       iEvent.getByLabel(gentopjet_sources[j], reco_gentopjets);
715 +
716 +       for (unsigned int i = 0; i < reco_gentopjets->size(); i++) {
717 +        
718 +         const reco::BasicJet  reco_gentopjet =  reco_gentopjets->at(i);
719 +         if(reco_gentopjet.pt() < gentopjet_ptmin) continue;
720 +         if(fabs(reco_gentopjet.eta()) > gentopjet_etamax) continue;
721 +
722 +         TopJet gentopjet;
723 +         gentopjet.set_charge(reco_gentopjet.charge());
724 +         gentopjet.set_pt(reco_gentopjet.pt());
725 +         gentopjet.set_eta(reco_gentopjet.eta());
726 +         gentopjet.set_phi(reco_gentopjet.phi());
727 +         gentopjet.set_energy(reco_gentopjet.energy());
728 +         gentopjet.set_numberOfDaughters(reco_gentopjet.numberOfDaughters());
729 +
730 +         for (unsigned int k = 0; k < reco_gentopjet.numberOfDaughters(); k++) {
731 +           Particle subjet_v4;
732 +           subjet_v4.set_pt(reco_gentopjet.daughter(k)->p4().pt());
733 +           subjet_v4.set_eta(reco_gentopjet.daughter(k)->p4().eta());
734 +           subjet_v4.set_phi(reco_gentopjet.daughter(k)->p4().phi());
735 +           subjet_v4.set_energy(reco_gentopjet.daughter(k)->p4().E());
736 +           gentopjet.add_subjet(subjet_v4);
737 +         }
738 +         gentopjets[j].push_back(gentopjet);
739 +       }
740 +     }
741 +   }
742 +
743     // ------------- photons -------------
744     if(doPhotons){
745       for(size_t j=0; j< photon_sources.size(); ++j){
# Line 496 | Line 752 | NtupleWriter::analyze(const edm::Event&
752         for (unsigned int i = 0; i < pat_photons.size(); ++i) {
753           pat::Photon pat_photon = pat_photons[i];
754           Photon ph;
755 <         ph.charge = 0;
756 <         ph.pt =  pat_photon.pt();
757 <         ph.eta =  pat_photon.eta();
758 <         ph.phi =  pat_photon.phi();
759 <         ph.energy =  pat_photon.energy();
760 <         ph.vertex_x = pat_photon.vertex().x();
761 <         ph.vertex_y = pat_photon.vertex().y();
762 <         ph.vertex_z = pat_photon.vertex().z();
763 <         ph.supercluster_eta = pat_photon.superCluster()->eta();
764 <         ph.supercluster_phi = pat_photon.superCluster()->phi();
765 < //       ph.neutralHadronIso = pat_photon.neutralHadronIso();
766 < //       ph.chargedHadronIso = pat_photon.chargedHadronIso();
767 <         ph.trackIso = pat_photon.trackIso();
755 >         ph.set_charge(0);
756 >         ph.set_pt( pat_photon.pt());
757 >         ph.set_eta( pat_photon.eta());
758 >         ph.set_phi( pat_photon.phi());
759 >         ph.set_energy( pat_photon.energy());
760 >         ph.set_vertex_x(pat_photon.vertex().x());
761 >         ph.set_vertex_y(pat_photon.vertex().y());
762 >         ph.set_vertex_z(pat_photon.vertex().z());
763 >         ph.set_supercluster_eta(pat_photon.superCluster()->eta());
764 >         ph.set_supercluster_phi(pat_photon.superCluster()->phi());
765 > //       ph.set_neutralHadronIso(pat_photon.neutralHadronIso());
766 > //       ph.set_chargedHadronIso(pat_photon.chargedHadronIso());
767 >         ph.set_trackIso(pat_photon.trackIso());
768           phs[j].push_back(ph);
769         }
770       }
# Line 528 | Line 784 | NtupleWriter::analyze(const edm::Event&
784         else{
785           pat::MET pat_met = pat_mets[0];
786          
787 <         met[j].pt=pat_met.pt();
788 <         met[j].phi=pat_met.phi();
789 <         met[j].mEtSig= pat_met.mEtSig();
787 >         met[j].set_pt(pat_met.pt());
788 >         met[j].set_phi(pat_met.phi());
789 >         met[j].set_mEtSig(pat_met.mEtSig());
790         }
791        
792       }
793     }
794  
795     // ------------- trigger -------------
796 <  
797 <   edm::InputTag triggerEvent = edm::InputTag("hltTriggerSummaryAOD");
798 <   edm::Handle< trigger::TriggerEvent > dummy_TriggerEvent;
799 <   iEvent.getByLabel( edm::InputTag(triggerEvent.label(), triggerEvent.instance()), dummy_TriggerEvent );
544 <  
545 <   const edm::Provenance *meta = dummy_TriggerEvent.provenance();
546 <   std::string nameProcess = meta->processName();
547 <   edm::InputTag triggerResultTag = edm::InputTag("TriggerResults");
548 <   triggerResultTag = edm::InputTag( triggerResultTag.label(), triggerResultTag.instance(), nameProcess );
549 <  
550 <   edm::Handle<edm::TriggerResults> trigger;
551 <   iEvent.getByLabel(triggerResultTag, trigger);
552 <   const edm::TriggerResults& trig = *(trigger.product());
553 <  
554 <   triggerResults.clear();
555 <   triggerNames.clear();
556 <   L1_prescale.clear();
557 <   HLT_prescale.clear();
558 <
559 <   edm::Service<edm::service::TriggerNamesService> tns;
560 <   std::vector<std::string> triggerNames_all;
561 <   tns->getTrigPaths(trig,triggerNames_all);
562 <
563 <   if (trig.size()!=triggerNames_all.size()) std::cout <<"ERROR: length of names and paths not the same: "<<triggerNames_all.size()<<","<<trig.size()<< std::endl;
564 <   for(unsigned int i=0; i<trig.size(); ++i){
565 <     std::vector<std::string>::const_iterator it = trigger_prefixes.begin();
566 <     for(; it!=trigger_prefixes.end(); ++it){
567 <       if(triggerNames_all[i].substr(0, it->size()) == *it)break;
568 <     }
569 <     if(it==trigger_prefixes.end()) continue;
570 <
571 <     //triggerResults.insert(std::pair<std::string, bool>(triggerNames[i],trig.accept(i)));
572 <     triggerResults.push_back(trig.accept(i));
573 <     if(newrun) triggerNames.push_back(triggerNames_all[i]);
574 <     if(isRealData){
575 <       std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
576 <       L1_prescale.push_back(pre.first);
577 <       HLT_prescale.push_back(pre.second);
578 <     }
579 <   }
580 < //    for(std::map<std::string, bool>::const_iterator iter = triggerResults.begin(); iter!=triggerResults.end(); iter++){
581 < //      std::cout << (*iter).first << "   " << (*iter).second << std::endl;
582 < //    }
583 <   newrun=false;
584 <
585 <   // ------------- generator info -------------
586 <  
587 <   if(doGenInfo){
588 <     genInfo.weights.clear();
589 <     genInfo.binningValues.clear();
590 <     genps.clear();
591 <
592 <     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
593 <     iEvent.getByLabel("generator", genEventInfoProduct);
594 <     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
595 <  
596 <     genInfo.binningValues = genEventInfo.binningValues();
597 <     genInfo.weights = genEventInfo.weights();
598 <     genInfo.alphaQCD = genEventInfo.alphaQCD();
599 <     genInfo.alphaQED = genEventInfo.alphaQED();
600 <     genInfo.qScale = genEventInfo.qScale();
796 >   if(doTrigger){
797 >     edm::InputTag triggerEvent = edm::InputTag("hltTriggerSummaryAOD");
798 >     edm::Handle< trigger::TriggerEvent > dummy_TriggerEvent;
799 >     iEvent.getByLabel( edm::InputTag(triggerEvent.label(), triggerEvent.instance()), dummy_TriggerEvent );
800      
801 <     const gen::PdfInfo* pdf = genEventInfo.pdf();
802 <     if(pdf){
803 <       genInfo.pdf_id1=pdf->id.first;
804 <       genInfo.pdf_id2=pdf->id.second;
805 <       genInfo.pdf_x1=pdf->x.first;
806 <       genInfo.pdf_x2=pdf->x.second;
807 <       genInfo.pdf_scalePDF=pdf->scalePDF;
808 <       genInfo.pdf_xPDF1=pdf->xPDF.first;
809 <       genInfo.pdf_xPDF2=pdf->xPDF.second;
810 <     }
811 <     else{
812 <       genInfo.pdf_id1=-999;
813 <       genInfo.pdf_id2=-999;
814 <       genInfo.pdf_x1=-999;
815 <       genInfo.pdf_x2=-999;
816 <       genInfo.pdf_scalePDF=-999;
817 <       genInfo.pdf_xPDF1=-999;
818 <       genInfo.pdf_xPDF2=-999;
819 <     }
820 <
821 <     edm::Handle<std::vector<PileupSummaryInfo> > pus;
822 <     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
823 <     genInfo.pileup_NumInteractions_intime=0;
824 <     genInfo.pileup_NumInteractions_ootbefore=0;
825 <     genInfo.pileup_NumInteractions_ootafter=0;
826 <     if(pus.isValid()){
827 <       genInfo.pileup_TrueNumInteractions = (float) pus->at(0).getTrueNumInteractions();
828 <       for(size_t i=0; i<pus->size(); ++i){
829 <         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
830 <            genInfo.pileup_NumInteractions_intime += pus->at(i).getPU_NumInteractions();
831 <         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
832 <            genInfo.pileup_NumInteractions_ootbefore += pus->at(i).getPU_NumInteractions();
833 <         }
834 <         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
835 <           genInfo.pileup_NumInteractions_ootafter += pus->at(i).getPU_NumInteractions();
836 <         }
837 <       }
838 <     }
839 <
840 <     edm::Handle<reco::GenParticleCollection> genPartColl;
642 <     iEvent.getByLabel(edm::InputTag("genParticles"), genPartColl);
643 <     int index=-1;
644 <     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
645 <       index++;
646 <      
647 <       //write out only top quarks and status 3 particles (works fine only for MadGraph)
648 <       if(abs(iter->pdgId())==6 || iter->status()==3){
649 <         GenParticle genp;
650 <         genp.charge = iter->charge();;
651 <         genp.pt = iter->p4().pt();
652 <         genp.eta = iter->p4().eta();
653 <         genp.phi = iter->p4().phi();
654 <         genp.energy = iter->p4().E();
655 <         genp.index =index;
656 <         genp.status = iter->status();
657 <         genp.pdgId = iter->pdgId();
658 <
659 <         genp.mother1=-1;
660 <         genp.mother2=-1;
661 <         genp.daughter1=-1;
662 <         genp.daughter2=-1;
663 <        
664 <         int nm=iter->numberOfMothers();
665 <         int nd=iter->numberOfDaughters();
666 <        
667 <         if (nm>0) genp.mother1 = iter->motherRef(0).key();
668 <         if (nm>1) genp.mother2 = iter->motherRef(nm-1).key();
669 <         if (nd>0) genp.daughter1 = iter->daughterRef(0).key();
670 <         if (nd>1) genp.daughter2 = iter->daughterRef(nd-1).key();
671 <        
672 <         genps.push_back(genp);
673 <       }
674 <     }
675 <
801 >     const edm::Provenance *meta = dummy_TriggerEvent.provenance();
802 >     std::string nameProcess = meta->processName();
803 >     edm::InputTag triggerResultTag = edm::InputTag("TriggerResults");
804 >     triggerResultTag = edm::InputTag( triggerResultTag.label(), triggerResultTag.instance(), nameProcess );
805 >    
806 >     edm::Handle<edm::TriggerResults> trigger;
807 >     iEvent.getByLabel(triggerResultTag, trigger);
808 >     const edm::TriggerResults& trig = *(trigger.product());
809 >    
810 >     triggerResults.clear();
811 >     triggerNames.clear();
812 > //      L1_prescale.clear();
813 > //      HLT_prescale.clear();
814 >    
815 >     edm::Service<edm::service::TriggerNamesService> tns;
816 >     std::vector<std::string> triggerNames_all;
817 >     tns->getTrigPaths(trig,triggerNames_all);
818 >    
819 >     if (trig.size()!=triggerNames_all.size()) std::cout <<"ERROR: length of names and paths not the same: "<<triggerNames_all.size()<<","<<trig.size()<< std::endl;
820 >     for(unsigned int i=0; i<trig.size(); ++i){
821 >       std::vector<std::string>::const_iterator it = trigger_prefixes.begin();
822 >       for(; it!=trigger_prefixes.end(); ++it){
823 >         if(triggerNames_all[i].substr(0, it->size()) == *it)break;
824 >       }
825 >       if(it==trigger_prefixes.end()) continue;
826 >      
827 >       //triggerResults.insert(std::pair<std::string, bool>(triggerNames[i],trig.accept(i)));
828 >       triggerResults.push_back(trig.accept(i));
829 >       if(newrun) triggerNames.push_back(triggerNames_all[i]);
830 > //        if(isRealData){
831 > //       std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
832 > //       L1_prescale.push_back(pre.first);
833 > //       HLT_prescale.push_back(pre.second);
834 > //       //std::cout <<  triggerNames_all[i] << " " << pre.first << " " <<pre.second << "   " << hlt_cfg.prescaleValue(iEvent, iSetup, triggerNames_all[i]) << std::endl;
835 > //        }
836 >     }
837 >     //    for(std::map<std::string, bool>::const_iterator iter = triggerResults.begin(); iter!=triggerResults.end(); iter++){
838 >     //      std::cout << (*iter).first << "   " << (*iter).second << std::endl;
839 >     //    }
840 >     newrun=false;
841     }
842  
843  
844     tr->Fill();
845 <
845 >   if(doLumiInfo)
846 >     previouslumiblockwasfilled=true;
847   }
848  
849  
# Line 685 | Line 851 | NtupleWriter::analyze(const edm::Event&
851   void
852   NtupleWriter::beginJob()
853   {
854 +  if(doLumiInfo){
855 +    totalRecLumi=0;
856 +    totalDelLumi=0;
857 +    previouslumiblockwasfilled=false;
858 +  }
859   }
860  
861   // ------------ method called once each job just after ending the event loop  ------------
# Line 700 | Line 871 | NtupleWriter::endJob()
871   void
872   NtupleWriter::beginRun(edm::Run const& iRun, edm::EventSetup const&  iSetup)
873   {
874 <  bool setup_changed = false;
875 <  hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
876 <  newrun=true;
877 <
878 <  edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
708 <  iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
709 <  JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
710 <  jecUnc = new JetCorrectionUncertainty(JetCorPar);
874 >  if(doTrigger){
875 >    //bool setup_changed = false;
876 >    //hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
877 >    newrun=true;
878 >  }
879  
880 +  if(doJets || doTopJets){
881 +    if(doJECUncertainty){
882 +      edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
883 +      iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
884 +      JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
885 +      jecUnc = new JetCorrectionUncertainty(JetCorPar);
886 +    }
887 +  }
888   }
889  
890   // ------------ method called when ending the processing of a run  ------------
891   void
892   NtupleWriter::endRun(edm::Run const&, edm::EventSetup const&)
893   {
894 +  if(doLumiInfo)
895 +    std::cout << "total integ. luminosity: " << totalDelLumi <<"(del) " << totalRecLumi << "(rec)" << std::endl;
896   }
897  
898   // ------------ method called when starting to processes a luminosity block  ------------
899   void
900 < NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&)
900 > NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const& lumi, edm::EventSetup const&)
901   {
902 +  if(doLumiInfo){
903 +    edm::Handle<LumiSummary> l;
904 +    lumi.getByLabel("lumiProducer", l);
905 +    
906 +    //add lumi of lumi blocks without any event to next lumiblock
907 +    if(previouslumiblockwasfilled){
908 +      intgRecLumi=0;
909 +      intgDelLumi=0;
910 +    }
911 +    previouslumiblockwasfilled=false;
912 +    
913 +    if (l.isValid()){;
914 +      intgRecLumi+=l->intgRecLumi()*6.37;
915 +      intgDelLumi+=l->intgDelLumi()*6.37;
916 +      totalRecLumi+=l->intgRecLumi()*6.37;
917 +      totalDelLumi+=l->intgDelLumi()*6.37;
918 +    }
919 +    //std::cout << "this lb: " <<l->intgRecLumi()*6.37 <<"   " << l->intgDelLumi()*6.37<<std::endl;
920 +    //std::cout << "summed:  "<< intgRecLumi << "   " << intgDelLumi << std::endl;
921 +  }
922   }
923  
924   // ------------ method called when ending the processing of a luminosity block  ------------

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines