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.9 by peiffer, Mon Apr 16 15:42:09 2012 UTC vs.
Revision 1.22 by peiffer, Wed Jul 25 09:56:57 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");
# Line 63 | 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);
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);
# Line 105 | 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 124 | Line 143 | NtupleWriter::NtupleWriter(const edm::Pa
143      }
144    }
145    if(doGenInfo){
146 +    genparticle_source= iConfig.getParameter<edm::InputTag>("genparticle_source");
147      tr->Branch("genInfo","GenInfo",&genInfo);
148      tr->Branch("GenParticles","std::vector<GenParticle>", &genps);
149    }
# Line 132 | Line 152 | NtupleWriter::NtupleWriter(const edm::Pa
152      //tr->Branch("triggerResults","std::map<std::string, bool>",&triggerResults);
153      tr->Branch("triggerNames", "std::vector<std::string>", &triggerNames);  
154      tr->Branch("triggerResults", "std::vector<bool>", &triggerResults);
155 <    tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
156 <    tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
155 >    //tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
156 >    //tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
157    }
158    newrun = true;
159   }
# Line 166 | Line 186 | NtupleWriter::analyze(const edm::Event&
186     luminosityBlock = iEvent.luminosityBlock();
187     isRealData      = iEvent.isRealData();
188  
189 <   if(isRealData){
190 <     edm::Handle<bool> bool_handle;
191 <     iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
192 <     HBHENoiseFilterResult = *(bool_handle.product());
193 <   }
194 <   else HBHENoiseFilterResult = false;
189 >   edm::Handle<double> m_rho;
190 >   iEvent.getByLabel(rho_source,m_rho);
191 >   rho=*m_rho;
192 >
193 > //    if(isRealData){
194 > //      edm::Handle<bool> bool_handle;
195 > //      iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
196 > //      HBHENoiseFilterResult = *(bool_handle.product());
197 > //    }
198 > //    else HBHENoiseFilterResult = false;
199  
200     // ------------- primary vertices and beamspot  -------------
201  
# Line 187 | Line 211 | NtupleWriter::analyze(const edm::Event&
211           reco::Vertex reco_pv = reco_pvs[i];
212  
213           PrimaryVertex pv;
214 <         pv.x =  reco_pv.x();
215 <         pv.y =  reco_pv.y();
216 <         pv.z =  reco_pv.z();
217 <         pv.nTracks =  reco_pv.nTracks();
218 <         //pv.isValid =  reco_pv.isValid();
219 <         pv.chi2 =  reco_pv.chi2();
220 <         pv.ndof =  reco_pv.ndof();      
214 >         pv.set_x( reco_pv.x());
215 >         pv.set_y( reco_pv.y());
216 >         pv.set_z( reco_pv.z());
217 >         pv.set_nTracks( reco_pv.nTracks());
218 >         //pv.set_isValid( reco_pv.isValid());
219 >         pv.set_chi2( reco_pv.chi2());
220 >         pv.set_ndof( reco_pv.ndof());  
221  
222           pvs[j].push_back(pv);
223         }
# Line 208 | Line 232 | NtupleWriter::analyze(const edm::Event&
232       beamspot_z0 = bsp.z0();
233     }
234  
235 + // ------------- generator info -------------
236 +  
237 +   if(doGenInfo){
238 +     genInfo.clear_weights();
239 +     genInfo.clear_binningValues();
240 +     genps.clear();
241 +
242 +     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
243 +     iEvent.getByLabel("generator", genEventInfoProduct);
244 +     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
245 +  
246 +     for(unsigned int k=0; k<genEventInfo.binningValues().size();++k){
247 +       genInfo.add_binningValue(genEventInfo.binningValues().at(k));
248 +     }
249 +     for(unsigned int k=0; k<genEventInfo.weights().size();++k){
250 +       genInfo.add_weight(genEventInfo.weights().at(k));
251 +     }
252 +     genInfo.set_alphaQCD(genEventInfo.alphaQCD());
253 +     genInfo.set_alphaQED(genEventInfo.alphaQED());
254 +     genInfo.set_qScale(genEventInfo.qScale());
255 +    
256 +     const gen::PdfInfo* pdf = genEventInfo.pdf();
257 +     if(pdf){
258 +       genInfo.set_pdf_id1(pdf->id.first);
259 +       genInfo.set_pdf_id2(pdf->id.second);
260 +       genInfo.set_pdf_x1(pdf->x.first);
261 +       genInfo.set_pdf_x2(pdf->x.second);
262 +       genInfo.set_pdf_scalePDF(pdf->scalePDF);
263 +       genInfo.set_pdf_xPDF1(pdf->xPDF.first);
264 +       genInfo.set_pdf_xPDF2(pdf->xPDF.second);
265 +     }
266 +     else{
267 +       genInfo.set_pdf_id1(-999);
268 +       genInfo.set_pdf_id2(-999);
269 +       genInfo.set_pdf_x1(-999);
270 +       genInfo.set_pdf_x2(-999);
271 +       genInfo.set_pdf_scalePDF(-999);
272 +       genInfo.set_pdf_xPDF1(-999);
273 +       genInfo.set_pdf_xPDF2(-999);
274 +     }
275 +
276 +     edm::Handle<std::vector<PileupSummaryInfo> > pus;
277 +     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
278 +     genInfo.set_pileup_NumInteractions_intime(0);
279 +     genInfo.set_pileup_NumInteractions_ootbefore(0);
280 +     genInfo.set_pileup_NumInteractions_ootafter(0);
281 +     if(pus.isValid()){
282 +       genInfo.set_pileup_TrueNumInteractions ( (float) pus->at(0).getTrueNumInteractions());
283 +       for(size_t i=0; i<pus->size(); ++i){
284 +         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
285 +           genInfo.set_pileup_NumInteractions_intime( genInfo.pileup_NumInteractions_intime() + pus->at(i).getPU_NumInteractions());
286 +         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
287 +           genInfo.set_pileup_NumInteractions_ootbefore( genInfo.pileup_NumInteractions_ootbefore() + pus->at(i).getPU_NumInteractions());
288 +         }
289 +         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
290 +           genInfo.set_pileup_NumInteractions_ootafter( genInfo.pileup_NumInteractions_ootafter() + pus->at(i).getPU_NumInteractions());
291 +         }
292 +       }
293 +     }
294 +
295 +     edm::Handle<reco::GenParticleCollection> genPartColl;
296 +     iEvent.getByLabel(genparticle_source, genPartColl);
297 +     int index=-1;
298 +     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
299 +       index++;
300 +      
301 +       //write out only top quarks and status 3 particles (works fine only for MadGraph)
302 +       if(abs(iter->pdgId())==6 || iter->status()==3 || doAllGenParticles){
303 +         GenParticle genp;
304 +         genp.set_charge(iter->charge());
305 +         genp.set_pt(iter->p4().pt());
306 +         genp.set_eta(iter->p4().eta());
307 +         genp.set_phi(iter->p4().phi());
308 +         genp.set_energy(iter->p4().E());
309 +         genp.set_index(index);
310 +         genp.set_status( iter->status());
311 +         genp.set_pdgId( iter->pdgId());
312 +
313 +         genp.set_mother1(-1);
314 +         genp.set_mother2(-1);
315 +         genp.set_daughter1(-1);
316 +         genp.set_daughter2(-1);
317 +        
318 +         int nm=iter->numberOfMothers();
319 +         int nd=iter->numberOfDaughters();
320 +
321 +        
322 +         if (nm>0) genp.set_mother1( iter->motherRef(0).key());
323 +         if (nm>1) genp.set_mother2( iter->motherRef(1).key());
324 +         if (nd>0) genp.set_daughter1( iter->daughterRef(0).key());
325 +         if (nd>1) genp.set_daughter2( iter->daughterRef(1).key());
326 +
327 +         genps.push_back(genp);
328 +       }
329 +     }
330 +   }
331 +
332     // ------------- electrons -------------  
333     if(doElectrons){
334 +
335 + //      edm::Handle<reco::ConversionCollection> hConversions;
336 + //      iEvent.getByLabel("allConversions", hConversions);
337 +
338 + //      edm::Handle<reco::BeamSpot> beamSpot;
339 + //      iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
340 + //      const reco::BeamSpot & bsp = *beamSpot;
341 +
342       for(size_t j=0; j< electron_sources.size(); ++j){
343         eles[j].clear();
344         edm::Handle< std::vector<pat::Electron> > ele_handle;
# Line 220 | Line 349 | NtupleWriter::analyze(const edm::Event&
349           pat::Electron pat_ele = pat_electrons[i];
350           Electron ele;
351          
352 <         ele.charge =  pat_ele.charge();
353 <         ele.pt =  pat_ele.pt();
354 <         ele.eta =  pat_ele.eta();
355 <         ele.phi =  pat_ele.phi();
356 <         ele.energy =  pat_ele.energy();
357 <         ele.vertex_x = pat_ele.vertex().x();
358 <         ele.vertex_y = pat_ele.vertex().y();
359 <         ele.vertex_z = pat_ele.vertex().z();
360 <         ele.supercluster_eta = pat_ele.superCluster()->eta();
361 <         ele.supercluster_phi = pat_ele.superCluster()->phi();
362 <         ele.dB = pat_ele.dB();
363 <         //ele.particleIso = pat_ele.particleIso();
364 <         ele.neutralHadronIso = pat_ele.neutralHadronIso();
365 <         ele.chargedHadronIso = pat_ele.chargedHadronIso();
366 <         ele.trackIso = pat_ele.trackIso();
367 <         ele.puChargedHadronIso = pat_ele.puChargedHadronIso();
368 <         ele.gsfTrack_trackerExpectedHitsInner_numberOfLostHits = pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits();
369 <         ele.gsfTrack_px= pat_ele.gsfTrack()->px();
370 <         ele.gsfTrack_py= pat_ele.gsfTrack()->py();
371 <         ele.gsfTrack_pz= pat_ele.gsfTrack()->pz();
372 <         ele.gsfTrack_vx= pat_ele.gsfTrack()->vx();
373 <         ele.gsfTrack_vy= pat_ele.gsfTrack()->vy();
374 <         ele.gsfTrack_vz= pat_ele.gsfTrack()->vz();
352 >         ele.set_charge( pat_ele.charge());
353 >         ele.set_pt( pat_ele.pt());
354 >         ele.set_eta( pat_ele.eta());
355 >         ele.set_phi( pat_ele.phi());
356 >         ele.set_energy( pat_ele.energy());
357 >         ele.set_vertex_x(pat_ele.vertex().x());
358 >         ele.set_vertex_y(pat_ele.vertex().y());
359 >         ele.set_vertex_z(pat_ele.vertex().z());
360 >         ele.set_supercluster_eta(pat_ele.superCluster()->eta());
361 >         ele.set_supercluster_phi(pat_ele.superCluster()->phi());
362 >         ele.set_dB(pat_ele.dB());
363 >         //ele.set_particleIso(pat_ele.particleIso());
364 >         ele.set_neutralHadronIso(pat_ele.neutralHadronIso());
365 >         ele.set_chargedHadronIso(pat_ele.chargedHadronIso());
366 >         ele.set_trackIso(pat_ele.trackIso());
367 >         ele.set_photonIso(pat_ele.photonIso());
368 >         ele.set_puChargedHadronIso(pat_ele.puChargedHadronIso());
369 >         ele.set_gsfTrack_trackerExpectedHitsInner_numberOfLostHits(pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits());
370 >         ele.set_gsfTrack_px( pat_ele.gsfTrack()->px());
371 >         ele.set_gsfTrack_py( pat_ele.gsfTrack()->py());
372 >         ele.set_gsfTrack_pz( pat_ele.gsfTrack()->pz());
373 >         ele.set_gsfTrack_vx( pat_ele.gsfTrack()->vx());
374 >         ele.set_gsfTrack_vy( pat_ele.gsfTrack()->vy());
375 >         ele.set_gsfTrack_vz( pat_ele.gsfTrack()->vz());
376 >         //ele.set_passconversionveto(!ConversionTools::hasMatchedConversion(pat_ele,hConversions,bsp.position()));
377 >         ele.set_passconversionveto(pat_ele.passConversionVeto());
378 >         ele.set_dEtaIn(pat_ele.deltaEtaSuperClusterTrackAtVtx());
379 >         ele.set_dPhiIn(pat_ele.deltaPhiSuperClusterTrackAtVtx());
380 >         ele.set_sigmaIEtaIEta(pat_ele.sigmaIetaIeta());
381 >         ele.set_HoverE(pat_ele.hadronicOverEm());
382 >         ele.set_fbrem(pat_ele.fbrem());
383 >         ele.set_EoverPIn(pat_ele.eSuperClusterOverP());
384 >         ele.set_EcalEnergy(pat_ele.ecalEnergy());
385 >         ele.set_mvaTrigV0(pat_ele.electronID("mvaTrigV0"));
386 >         ele.set_mvaNonTrigV0(pat_ele.electronID("mvaNonTrigV0"));
387 >         float AEff03 = 0.00;
388 >         if(isRealData){
389 >           AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAData2011);
390 >         }else{
391 >           AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAFall11MC);
392 >         }
393 >         ele.set_AEff(AEff03);
394 >
395           eles[j].push_back(ele);
396         }
397       }
# Line 261 | Line 410 | NtupleWriter::analyze(const edm::Event&
410           pat::Muon pat_mu = pat_muons[i];
411  
412           Muon mu;
413 <         mu.charge =  pat_mu.charge();
414 <         mu.pt =  pat_mu.pt();
415 <         mu.eta =  pat_mu.eta();
416 <         mu.phi =  pat_mu.phi();
417 <         mu.energy =  pat_mu.energy();
418 <         mu.vertex_x = pat_mu.vertex().x();
419 <         mu.vertex_y = pat_mu.vertex().y();
420 <         mu.vertex_z = pat_mu.vertex().z();
421 <         mu.dB = pat_mu.dB();
422 <         //mu.particleIso = pat_mu.particleIso();
423 <         mu.neutralHadronIso = pat_mu.neutralHadronIso();
424 <         mu.chargedHadronIso = pat_mu.chargedHadronIso();
425 <         mu.trackIso = pat_mu.trackIso();
426 <         mu.puChargedHadronIso = pat_mu.puChargedHadronIso();
427 <         mu.isGlobalMuon = pat_mu.isGlobalMuon();
428 <         mu.isStandAloneMuon = pat_mu.isStandAloneMuon();
429 <         mu.isTrackerMuon = pat_mu.isTrackerMuon();
430 <         mu.numberOfMatchedStations = pat_mu.numberOfMatchedStations();
413 >         mu.set_charge( pat_mu.charge());
414 >         mu.set_pt( pat_mu.pt());
415 >         mu.set_eta( pat_mu.eta());
416 >         mu.set_phi( pat_mu.phi());
417 >         mu.set_energy( pat_mu.energy());
418 >         mu.set_vertex_x ( pat_mu.vertex().x());
419 >         mu.set_vertex_y ( pat_mu.vertex().y());
420 >         mu.set_vertex_z ( pat_mu.vertex().z());
421 >         mu.set_dB ( pat_mu.dB());
422 >         //mu.particleIso ( pat_mu.particleIso());
423 >         mu.set_neutralHadronIso ( pat_mu.neutralHadronIso());
424 >         mu.set_chargedHadronIso ( pat_mu.chargedHadronIso());
425 >         mu.set_trackIso ( pat_mu.trackIso());
426 >         mu.set_photonIso ( pat_mu.photonIso());
427 >         mu.set_puChargedHadronIso ( pat_mu.puChargedHadronIso());
428 >         mu.set_isGlobalMuon ( pat_mu.isGlobalMuon());
429 >         mu.set_isPFMuon ( pat_mu.isPFMuon());
430 >         mu.set_isStandAloneMuon ( pat_mu.isStandAloneMuon());
431 >         mu.set_isTrackerMuon ( pat_mu.isTrackerMuon());
432 >         mu.set_numberOfMatchedStations ( pat_mu.numberOfMatchedStations());
433           reco::TrackRef globalTrack = pat_mu.globalTrack();
434           if(!globalTrack.isNull()){
435 <           mu.globalTrack_chi2 = globalTrack->chi2();
436 <           mu.globalTrack_ndof = globalTrack->ndof();
437 <           mu.globalTrack_d0 = globalTrack->d0();        
438 <           mu.globalTrack_d0Error = globalTrack->d0Error();
439 <           mu.globalTrack_numberOfValidHits = globalTrack->numberOfValidHits();
440 <           mu.globalTrack_numberOfLostHits = globalTrack->numberOfLostHits();
435 >           mu.set_globalTrack_chi2 ( globalTrack->chi2());
436 >           mu.set_globalTrack_ndof ( globalTrack->ndof());
437 >           mu.set_globalTrack_d0 ( globalTrack->d0());  
438 >           mu.set_globalTrack_d0Error ( globalTrack->d0Error());
439 >           mu.set_globalTrack_numberOfValidHits ( globalTrack->numberOfValidHits());
440 >           mu.set_globalTrack_numberOfLostHits ( globalTrack->numberOfLostHits());
441 >           mu.set_globalTrack_numberOfValidMuonHits(globalTrack->hitPattern().numberOfValidMuonHits() );
442           }
443           else{
444 <           mu.globalTrack_chi2 = 0;
445 <           mu.globalTrack_ndof = 0;
446 <           mu.globalTrack_d0 = 0;
447 <           mu.globalTrack_d0Error = 0;
448 <           mu.globalTrack_numberOfValidHits = 0;
449 <           mu.globalTrack_numberOfLostHits = 0;
444 >           mu.set_globalTrack_chi2 ( 0);
445 >           mu.set_globalTrack_ndof ( 0);
446 >           mu.set_globalTrack_d0 ( 0);
447 >           mu.set_globalTrack_d0Error ( 0);
448 >           mu.set_globalTrack_numberOfValidHits ( 0);
449 >           mu.set_globalTrack_numberOfLostHits ( 0);
450           }
451           reco::TrackRef innerTrack = pat_mu.innerTrack();
452           if(!innerTrack.isNull()){
453 <           mu.innerTrack_chi2 = innerTrack->chi2();
454 <           mu.innerTrack_ndof = innerTrack->ndof();
455 <           mu.innerTrack_d0 = innerTrack->d0();  
456 <           mu.innerTrack_d0Error = innerTrack->d0Error();
457 <           mu.innerTrack_numberOfValidHits = innerTrack->numberOfValidHits();
458 <           mu.innerTrack_numberOfLostHits = innerTrack->numberOfLostHits();
453 >           mu.set_innerTrack_chi2 ( innerTrack->chi2());
454 >           mu.set_innerTrack_ndof ( innerTrack->ndof());
455 >           mu.set_innerTrack_d0 ( innerTrack->d0());    
456 >           mu.set_innerTrack_d0Error ( innerTrack->d0Error());
457 >           mu.set_innerTrack_numberOfValidHits ( innerTrack->numberOfValidHits());
458 >           mu.set_innerTrack_numberOfLostHits ( innerTrack->numberOfLostHits());
459 >           mu.set_innerTrack_trackerLayersWithMeasurement ( innerTrack->hitPattern().trackerLayersWithMeasurement());
460 >           mu.set_innerTrack_numberOfValidPixelHits ( innerTrack->hitPattern().numberOfValidPixelHits());
461           }
462           else{
463 <           mu.innerTrack_chi2 = 0;
464 <           mu.innerTrack_ndof = 0;
465 <           mu.innerTrack_d0 = 0;
466 <           mu.innerTrack_d0Error = 0;
467 <           mu.innerTrack_numberOfValidHits = 0;
468 <           mu.innerTrack_numberOfLostHits = 0;
463 >           mu.set_innerTrack_chi2 ( 0);
464 >           mu.set_innerTrack_ndof ( 0);
465 >           mu.set_innerTrack_d0 ( 0);
466 >           mu.set_innerTrack_d0Error ( 0);
467 >           mu.set_innerTrack_numberOfValidHits ( 0);
468 >           mu.set_innerTrack_numberOfLostHits ( 0);
469 >           mu.set_innerTrack_trackerLayersWithMeasurement ( 0);
470 >           mu.set_innerTrack_numberOfValidPixelHits ( 0);
471           }
472           reco::TrackRef outerTrack = pat_mu.outerTrack();
473           if(!outerTrack.isNull()){
474 <           mu.outerTrack_chi2 = outerTrack->chi2();
475 <           mu.outerTrack_ndof = outerTrack->ndof();
476 <           mu.outerTrack_d0 = outerTrack->d0();  
477 <           mu.outerTrack_d0Error = outerTrack->d0Error();
478 <           mu.outerTrack_numberOfValidHits = outerTrack->numberOfValidHits();
479 <           mu.outerTrack_numberOfLostHits = outerTrack->numberOfLostHits();
474 >           mu.set_outerTrack_chi2 ( outerTrack->chi2());
475 >           mu.set_outerTrack_ndof ( outerTrack->ndof());
476 >           mu.set_outerTrack_d0 ( outerTrack->d0());    
477 >           mu.set_outerTrack_d0Error ( outerTrack->d0Error());
478 >           mu.set_outerTrack_numberOfValidHits ( outerTrack->numberOfValidHits());
479 >           mu.set_outerTrack_numberOfLostHits ( outerTrack->numberOfLostHits());
480           }
481           else{
482 <           mu.outerTrack_chi2 = 0;
483 <           mu.outerTrack_ndof = 0;
484 <           mu.outerTrack_d0 = 0;
485 <           mu.outerTrack_d0Error = 0;
486 <           mu.outerTrack_numberOfValidHits = 0;
487 <           mu.outerTrack_numberOfLostHits = 0;
482 >           mu.set_outerTrack_chi2 ( 0);
483 >           mu.set_outerTrack_ndof ( 0);
484 >           mu.set_outerTrack_d0 ( 0);
485 >           mu.set_outerTrack_d0Error ( 0);
486 >           mu.set_outerTrack_numberOfValidHits ( 0);
487 >           mu.set_outerTrack_numberOfLostHits ( 0);
488           }
489  
490           mus[j].push_back(mu);
# Line 352 | Line 508 | NtupleWriter::analyze(const edm::Event&
508           if(fabs(pat_tau.eta()) > tau_etamax) continue;
509  
510           Tau tau;
511 <         tau.charge =  pat_tau.charge();
512 <         tau.pt =  pat_tau.pt();
513 <         tau.eta =  pat_tau.eta();
514 <         tau.phi =  pat_tau.phi();
515 <         tau.energy =  pat_tau.energy();
516 <         tau.decayModeFinding = pat_tau.tauID("decayModeFinding")>0.5;
517 <         tau.byVLooseCombinedIsolationDeltaBetaCorr  = pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5;
518 <         tau.byLooseCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5;
519 <         tau.byMediumCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5;
520 <         tau.byTightCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5;
521 <         tau.againstElectronLoose  = pat_tau.tauID("againstElectronLoose")>0.5;
522 <         tau.againstElectronMedium = pat_tau.tauID("againstElectronMedium")>0.5;
523 <         tau.againstElectronTight = pat_tau.tauID("againstElectronTight")>0.5;
524 <         tau.againstElectronMVA  = pat_tau.tauID("againstElectronMVA")>0.5;
525 <         tau.againstMuonLoose = pat_tau.tauID("againstMuonLoose")>0.5;
526 <         tau.againstMuonMedium = pat_tau.tauID("againstMuonMedium")>0.5;
527 <         tau.againstMuonTight = pat_tau.tauID("againstMuonTight")>0.5;
528 <        
511 >         tau.set_charge( pat_tau.charge());
512 >         tau.set_pt( pat_tau.pt());
513 >         tau.set_eta( pat_tau.eta());
514 >         tau.set_phi( pat_tau.phi());
515 >         tau.set_energy( pat_tau.energy());
516 >         tau.set_decayModeFinding ( pat_tau.tauID("decayModeFinding")>0.5);
517 >         tau.set_byVLooseCombinedIsolationDeltaBetaCorr  ( pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5);
518 >         tau.set_byLooseCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5);
519 >         tau.set_byMediumCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5);
520 >         tau.set_byTightCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5);
521 >         tau.set_againstElectronLoose  ( pat_tau.tauID("againstElectronLoose")>0.5);
522 >         tau.set_againstElectronMedium ( pat_tau.tauID("againstElectronMedium")>0.5);
523 >         tau.set_againstElectronTight ( pat_tau.tauID("againstElectronTight")>0.5);
524 >         tau.set_againstElectronMVA  ( pat_tau.tauID("againstElectronMVA")>0.5);
525 >         tau.set_againstMuonLoose ( pat_tau.tauID("againstMuonLoose")>0.5);
526 >         tau.set_againstMuonMedium ( pat_tau.tauID("againstMuonMedium")>0.5);
527 >         tau.set_againstMuonTight ( pat_tau.tauID("againstMuonTight")>0.5);
528 >
529 > //       reco::PFCandidateRef leadPFCand = pat_tau.leadPFCand();
530 > //       if(!leadPFCand.isNull()){
531 > //         tau.set_leadPFCand_px ( leadPFCand->px());
532 > //         tau.set_leadPFCand_py ( leadPFCand->py());
533 > //         tau.set_leadPFCand_pz ( leadPFCand->pz());
534 > //       }
535 > //       else{
536 > //         tau.set_leadPFCand_px ( 0);
537 > //         tau.set_leadPFCand_py ( 0);
538 > //         tau.set_leadPFCand_pz ( 0);
539 > //       }
540           taus[j].push_back(tau);
541         }
542       }
# Line 396 | Line 563 | NtupleWriter::analyze(const edm::Event&
563   //       }
564  
565           Jet jet;
566 <         jet.charge = pat_jet.charge();
567 <         jet.pt = pat_jet.pt();
568 <         jet.eta = pat_jet.eta();
569 <         jet.phi = pat_jet.phi();
570 <         jet.energy = pat_jet.energy();
571 <         jet.numberOfDaughters =pat_jet.numberOfDaughters();
566 >         jet.set_charge(pat_jet.charge());
567 >         jet.set_pt(pat_jet.pt());
568 >         jet.set_eta(pat_jet.eta());
569 >         jet.set_phi(pat_jet.phi());
570 >         jet.set_energy(pat_jet.energy());
571 >         jet.set_numberOfDaughters (pat_jet.numberOfDaughters());
572           const reco::TrackRefVector&  jettracks = pat_jet.associatedTracks();
573 <         jet.nTracks = jettracks.size();
574 <         jet.jetArea = pat_jet.jetArea();
575 <         jet.pileup = pat_jet.pileup();
573 >         jet.set_nTracks ( jettracks.size());
574 >         jet.set_jetArea(pat_jet.jetArea());
575 >         jet.set_pileup(pat_jet.pileup());
576           if(pat_jet.isPFJet()){
577 <           jet.neutralEmEnergyFraction =pat_jet.neutralEmEnergyFraction();
578 <           jet.neutralHadronEnergyFraction =pat_jet.neutralHadronEnergyFraction();
579 <           jet.chargedEmEnergyFraction =pat_jet.chargedEmEnergyFraction();
580 <           jet.chargedHadronEnergyFraction =pat_jet.chargedHadronEnergyFraction();
581 <           jet.muonEnergyFraction =pat_jet.muonEnergyFraction();
582 <           jet.photonEnergyFraction =pat_jet.photonEnergyFraction();
583 <           jet.chargedMultiplicity =pat_jet.chargedMultiplicity();
584 <           jet.neutralMultiplicity =pat_jet.neutralMultiplicity();
585 <           jet.muonMultiplicity =pat_jet.muonMultiplicity();
586 <           jet.electronMultiplicity =pat_jet.electronMultiplicity();
587 <           jet.photonMultiplicity =pat_jet.photonMultiplicity();
588 <         }
589 <
590 <         jecUnc->setJetEta(pat_jet.eta());
591 <         jecUnc->setJetPt(pat_jet.pt());
592 <         jet.JEC_uncertainty = jecUnc->getUncertainty(true);
593 <
594 <         jet.btag_simpleSecondaryVertexHighEff=pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
595 <         jet.btag_simpleSecondaryVertexHighPur=pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
596 <         jet.btag_combinedSecondaryVertex=pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags");
597 <         jet.btag_combinedSecondaryVertexMVA=pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
598 <         jet.btag_jetBProbability=pat_jet.bDiscriminator("jetBProbabilityBJetTags");
599 <         jet.btag_jetProbability=pat_jet.bDiscriminator("jetProbabilityBJetTags");
577 >           jet.set_neutralEmEnergyFraction (pat_jet.neutralEmEnergyFraction());
578 >           jet.set_neutralHadronEnergyFraction (pat_jet.neutralHadronEnergyFraction());
579 >           jet.set_chargedEmEnergyFraction (pat_jet.chargedEmEnergyFraction());
580 >           jet.set_chargedHadronEnergyFraction (pat_jet.chargedHadronEnergyFraction());
581 >           jet.set_muonEnergyFraction (pat_jet.muonEnergyFraction());
582 >           jet.set_photonEnergyFraction (pat_jet.photonEnergyFraction());
583 >           jet.set_chargedMultiplicity (pat_jet.chargedMultiplicity());
584 >           jet.set_neutralMultiplicity (pat_jet.neutralMultiplicity());
585 >           jet.set_muonMultiplicity (pat_jet.muonMultiplicity());
586 >           jet.set_electronMultiplicity (pat_jet.electronMultiplicity());
587 >           jet.set_photonMultiplicity (pat_jet.photonMultiplicity());
588 >         }
589 >         if(doJECUncertainty){
590 >           jecUnc->setJetEta(pat_jet.eta());
591 >           jecUnc->setJetPt(pat_jet.pt());
592 >           jet.set_JEC_uncertainty(jecUnc->getUncertainty(true));
593 >         }
594 >         jet.set_JEC_factor_raw(pat_jet.jecFactor("Uncorrected"));
595 >        
596 >         jet.set_btag_simpleSecondaryVertexHighEff(pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
597 >         jet.set_btag_simpleSecondaryVertexHighPur(pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
598 >         jet.set_btag_combinedSecondaryVertex(pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags"));
599 >         jet.set_btag_combinedSecondaryVertexMVA(pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
600 >         jet.set_btag_jetBProbability(pat_jet.bDiscriminator("jetBProbabilityBJetTags"));
601 >         jet.set_btag_jetProbability(pat_jet.bDiscriminator("jetProbabilityBJetTags"));
602  
603 +        
604 +         const reco::GenJet *genj = pat_jet.genJet();
605 +         if(genj){
606 +           jet.set_genjet_pt(genj->pt());
607 +           jet.set_genjet_eta(genj->eta());
608 +           jet.set_genjet_phi(genj->phi());
609 +           jet.set_genjet_energy(genj->energy());
610 +           if(doAllGenParticles){
611 +             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
612 +             for(unsigned int l = 0; l<jetgenps.size(); ++l){
613 +               for(unsigned int k=0; k< genps.size(); ++k){
614 +                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
615 +                   jet.add_genparticles_index(genps[k].index());
616 +                   break;
617 +                 }
618 +               }
619 +             }
620 +             if(jet.genparticles_indices().size()!= jetgenps.size())
621 +               std::cout << "WARNING: Found only " << jet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this jet"<<std::endl;
622 +           }
623 +          
624 +         }
625 +        
626           jets[j].push_back(jet);
627         }
628       }
# Line 447 | Line 639 | NtupleWriter::analyze(const edm::Event&
639         iEvent.getByLabel(topjet_sources[j], pat_topjets);
640  
641         for (unsigned int i = 0; i < pat_topjets->size(); i++) {
642 +
643           const pat::Jet  pat_topjet =  * dynamic_cast<pat::Jet const *>(&pat_topjets->at(i));
644           if(pat_topjet.pt() < topjet_ptmin) continue;
645           if(fabs(pat_topjet.eta()) > topjet_etamax) continue;
646  
647           TopJet topjet;
648 <         topjet.charge = pat_topjet.charge();
649 <         topjet.pt = pat_topjet.pt();
650 <         topjet.eta = pat_topjet.eta();
651 <         topjet.phi = pat_topjet.phi();
652 <         topjet.energy = pat_topjet.energy();
653 <         topjet.numberOfDaughters =pat_topjet.numberOfDaughters();
648 >         topjet.set_charge(pat_topjet.charge());
649 >         topjet.set_pt(pat_topjet.pt());
650 >         topjet.set_eta(pat_topjet.eta());
651 >         topjet.set_phi(pat_topjet.phi());
652 >         topjet.set_energy(pat_topjet.energy());
653 >         topjet.set_numberOfDaughters(pat_topjet.numberOfDaughters());
654           const reco::TrackRefVector&  topjettracks = pat_topjet.associatedTracks();
655 <         topjet.nTracks = topjettracks.size();
656 <         topjet.jetArea = pat_topjet.jetArea();
657 <         topjet.pileup = pat_topjet.pileup();
658 < //       topjet.neutralEmEnergyFraction =pat_topjet.neutralEmEnergyFraction();
659 < //       topjet.neutralHadronEnergyFraction =pat_topjet.neutralHadronEnergyFraction();
660 < //       topjet.chargedEmEnergyFraction =pat_topjet.chargedEmEnergyFraction();
661 < //       topjet.chargedHadronEnergyFraction =pat_topjet.chargedHadronEnergyFraction();
662 < //       topjet.muonEnergyFraction =pat_topjet.muonEnergyFraction();
663 < //       topjet.photonEnergyFraction =pat_topjet.photonEnergyFraction();
664 < //       topjet.chargedMultiplicity =pat_topjet.chargedMultiplicity();
665 < //       topjet.neutralMultiplicity =pat_topjet.neutralMultiplicity();
666 < //       topjet.muonMultiplicity =pat_topjet.muonMultiplicity();
667 < //       topjet.electronMultiplicity =pat_topjet.electronMultiplicity();
668 < //       topjet.photonMultiplicity =pat_topjet.photonMultiplicity();
669 <
670 <         jecUnc->setJetEta(pat_topjet.eta());
671 <         jecUnc->setJetPt(pat_topjet.pt());
672 <         topjet.JEC_uncertainty = jecUnc->getUncertainty(true);
673 <
674 <         topjet.btag_simpleSecondaryVertexHighEff=pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
675 <         topjet.btag_simpleSecondaryVertexHighPur=pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
676 <         topjet.btag_combinedSecondaryVertex=pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags");
677 <         topjet.btag_combinedSecondaryVertexMVA=pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
678 <         topjet.btag_jetBProbability=pat_topjet.bDiscriminator("jetBProbabilityBJetTags");
679 <         topjet.btag_jetProbability=pat_topjet.bDiscriminator("jetProbabilityBJetTags");
655 >         topjet.set_nTracks( topjettracks.size());
656 >         topjet.set_jetArea( pat_topjet.jetArea());
657 >         topjet.set_pileup( pat_topjet.pileup());
658 > //       topjet.set_neutralEmEnergyFraction(pat_topjet.neutralEmEnergyFraction());
659 > //       topjet.set_neutralHadronEnergyFraction(pat_topjet.neutralHadronEnergyFraction());
660 > //       topjet.set_chargedEmEnergyFraction(pat_topjet.chargedEmEnergyFraction());
661 > //       topjet.set_chargedHadronEnergyFraction(pat_topjet.chargedHadronEnergyFraction());
662 > //       topjet.set_muonEnergyFraction(pat_topjet.muonEnergyFraction());
663 > //       topjet.set_photonEnergyFraction(pat_topjet.photonEnergyFraction());
664 > //       topjet.set_chargedMultiplicity(pat_topjet.chargedMultiplicity());
665 > //       topjet.set_neutralMultiplicity(pat_topjet.neutralMultiplicity());
666 > //       topjet.set_muonMultiplicity(pat_topjet.muonMultiplicity());
667 > //       topjet.set_electronMultiplicity(pat_topjet.electronMultiplicity());
668 > //       topjet.set_photonMultiplicity(pat_topjet.photonMultiplicity());
669 >         if(doJECUncertainty){
670 >           jecUnc->setJetEta(pat_topjet.eta());
671 >           jecUnc->setJetPt(pat_topjet.pt());
672 >           topjet.set_JEC_uncertainty( jecUnc->getUncertainty(true));
673 >         }
674 >         topjet.set_JEC_factor_raw( pat_topjet.jecFactor("Uncorrected"));
675 >
676 >         topjet.set_btag_simpleSecondaryVertexHighEff(pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
677 >         topjet.set_btag_simpleSecondaryVertexHighPur(pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
678 >         topjet.set_btag_combinedSecondaryVertex(pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags"));
679 >         topjet.set_btag_combinedSecondaryVertexMVA(pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
680 >         topjet.set_btag_jetBProbability(pat_topjet.bDiscriminator("jetBProbabilityBJetTags"));
681 >         topjet.set_btag_jetProbability(pat_topjet.bDiscriminator("jetProbabilityBJetTags"));
682 >
683 >         /*
684 >         const reco::GenJet *genj = pat_topjet.genJet();
685 >         if(genj){
686 >           topjet.set_genjet_pt ( genj->pt());
687 >           topjet.set_genjet_eta ( genj->eta());
688 >           topjet.set_genjet_phi ( genj->phi());
689 >           topjet.set_genjet_energy ( genj->energy());
690 >           if(doAllGenParticles){
691 >             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
692 >             for(unsigned int l = 0; l<jetgenps.size(); ++l){
693 >               for(unsigned int k=0; k< genps.size(); ++k){
694 >                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
695 >                   topjet.add_genparticles_index(genps[k].index());
696 >                   break;
697 >                 }
698 >               }
699 >             }
700 >             if(topjet.genparticles_indices().size()!= jetgenps.size())
701 >               std::cout << "WARNING: Found only " << topjet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this topjet"<<std::endl;
702 >           }
703 >         }
704 >         */
705  
706           for (unsigned int k = 0; k < pat_topjet.numberOfDaughters(); k++) {
707             Particle subjet_v4;
708 <           subjet_v4.pt = pat_topjet.daughter(k)->p4().pt();
709 <           subjet_v4.eta = pat_topjet.daughter(k)->p4().eta();
710 <           subjet_v4.phi = pat_topjet.daughter(k)->p4().phi();
711 <           subjet_v4.energy = pat_topjet.daughter(k)->p4().E();
712 <           topjet.subjets.push_back(subjet_v4);
708 >           subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt());
709 >           subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta());
710 >           subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi());
711 >           subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E());
712 >           topjet.add_subjet(subjet_v4);
713           }
714           topjets[j].push_back(topjet);
715         }
716       }
717     }
718  
719 +
720 +   // ------------- generator top jets -------------
721 +   if(doGenTopJets){
722 +     for(size_t j=0; j< gentopjet_sources.size(); ++j){
723 +      
724 +       gentopjets[j].clear();
725 +
726 +       edm::Handle<reco::BasicJetCollection> reco_gentopjets;
727 +       //edm::Handle<std::vector<reco::Jet> > reco_gentopjets;
728 +       iEvent.getByLabel(gentopjet_sources[j], reco_gentopjets);
729 +
730 +       for (unsigned int i = 0; i < reco_gentopjets->size(); i++) {
731 +        
732 +         const reco::BasicJet  reco_gentopjet =  reco_gentopjets->at(i);
733 +         if(reco_gentopjet.pt() < gentopjet_ptmin) continue;
734 +         if(fabs(reco_gentopjet.eta()) > gentopjet_etamax) continue;
735 +
736 +         TopJet gentopjet;
737 +         gentopjet.set_charge(reco_gentopjet.charge());
738 +         gentopjet.set_pt(reco_gentopjet.pt());
739 +         gentopjet.set_eta(reco_gentopjet.eta());
740 +         gentopjet.set_phi(reco_gentopjet.phi());
741 +         gentopjet.set_energy(reco_gentopjet.energy());
742 +         gentopjet.set_numberOfDaughters(reco_gentopjet.numberOfDaughters());
743 +
744 +         for (unsigned int k = 0; k < reco_gentopjet.numberOfDaughters(); k++) {
745 +           Particle subjet_v4;
746 +           subjet_v4.set_pt(reco_gentopjet.daughter(k)->p4().pt());
747 +           subjet_v4.set_eta(reco_gentopjet.daughter(k)->p4().eta());
748 +           subjet_v4.set_phi(reco_gentopjet.daughter(k)->p4().phi());
749 +           subjet_v4.set_energy(reco_gentopjet.daughter(k)->p4().E());
750 +           gentopjet.add_subjet(subjet_v4);
751 +         }
752 +         gentopjets[j].push_back(gentopjet);
753 +       }
754 +     }
755 +   }
756 +
757     // ------------- photons -------------
758     if(doPhotons){
759       for(size_t j=0; j< photon_sources.size(); ++j){
# Line 510 | Line 766 | NtupleWriter::analyze(const edm::Event&
766         for (unsigned int i = 0; i < pat_photons.size(); ++i) {
767           pat::Photon pat_photon = pat_photons[i];
768           Photon ph;
769 <         ph.charge = 0;
770 <         ph.pt =  pat_photon.pt();
771 <         ph.eta =  pat_photon.eta();
772 <         ph.phi =  pat_photon.phi();
773 <         ph.energy =  pat_photon.energy();
774 <         ph.vertex_x = pat_photon.vertex().x();
775 <         ph.vertex_y = pat_photon.vertex().y();
776 <         ph.vertex_z = pat_photon.vertex().z();
777 <         ph.supercluster_eta = pat_photon.superCluster()->eta();
778 <         ph.supercluster_phi = pat_photon.superCluster()->phi();
779 < //       ph.neutralHadronIso = pat_photon.neutralHadronIso();
780 < //       ph.chargedHadronIso = pat_photon.chargedHadronIso();
781 <         ph.trackIso = pat_photon.trackIso();
769 >         ph.set_charge(0);
770 >         ph.set_pt( pat_photon.pt());
771 >         ph.set_eta( pat_photon.eta());
772 >         ph.set_phi( pat_photon.phi());
773 >         ph.set_energy( pat_photon.energy());
774 >         ph.set_vertex_x(pat_photon.vertex().x());
775 >         ph.set_vertex_y(pat_photon.vertex().y());
776 >         ph.set_vertex_z(pat_photon.vertex().z());
777 >         ph.set_supercluster_eta(pat_photon.superCluster()->eta());
778 >         ph.set_supercluster_phi(pat_photon.superCluster()->phi());
779 > //       ph.set_neutralHadronIso(pat_photon.neutralHadronIso());
780 > //       ph.set_chargedHadronIso(pat_photon.chargedHadronIso());
781 >         ph.set_trackIso(pat_photon.trackIso());
782           phs[j].push_back(ph);
783         }
784       }
# Line 542 | Line 798 | NtupleWriter::analyze(const edm::Event&
798         else{
799           pat::MET pat_met = pat_mets[0];
800          
801 <         met[j].pt=pat_met.pt();
802 <         met[j].phi=pat_met.phi();
803 <         met[j].mEtSig= pat_met.mEtSig();
801 >         met[j].set_pt(pat_met.pt());
802 >         met[j].set_phi(pat_met.phi());
803 >         met[j].set_mEtSig(pat_met.mEtSig());
804         }
805        
806       }
# Line 567 | Line 823 | NtupleWriter::analyze(const edm::Event&
823      
824       triggerResults.clear();
825       triggerNames.clear();
826 <     L1_prescale.clear();
827 <     HLT_prescale.clear();
826 > //      L1_prescale.clear();
827 > //      HLT_prescale.clear();
828      
829       edm::Service<edm::service::TriggerNamesService> tns;
830       std::vector<std::string> triggerNames_all;
# Line 585 | Line 841 | NtupleWriter::analyze(const edm::Event&
841         //triggerResults.insert(std::pair<std::string, bool>(triggerNames[i],trig.accept(i)));
842         triggerResults.push_back(trig.accept(i));
843         if(newrun) triggerNames.push_back(triggerNames_all[i]);
844 <       if(isRealData){
845 <         std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
846 <         L1_prescale.push_back(pre.first);
847 <         HLT_prescale.push_back(pre.second);
848 <       }
844 > //        if(isRealData){
845 > //       std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
846 > //       L1_prescale.push_back(pre.first);
847 > //       HLT_prescale.push_back(pre.second);
848 > //       //std::cout <<  triggerNames_all[i] << " " << pre.first << " " <<pre.second << "   " << hlt_cfg.prescaleValue(iEvent, iSetup, triggerNames_all[i]) << std::endl;
849 > //        }
850       }
851       //    for(std::map<std::string, bool>::const_iterator iter = triggerResults.begin(); iter!=triggerResults.end(); iter++){
852       //      std::cout << (*iter).first << "   " << (*iter).second << std::endl;
# Line 597 | Line 854 | NtupleWriter::analyze(const edm::Event&
854       newrun=false;
855     }
856  
600   // ------------- generator info -------------
601  
602   if(doGenInfo){
603     genInfo.weights.clear();
604     genInfo.binningValues.clear();
605     genps.clear();
606
607     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
608     iEvent.getByLabel("generator", genEventInfoProduct);
609     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
610  
611     genInfo.binningValues = genEventInfo.binningValues();
612     genInfo.weights = genEventInfo.weights();
613     genInfo.alphaQCD = genEventInfo.alphaQCD();
614     genInfo.alphaQED = genEventInfo.alphaQED();
615     genInfo.qScale = genEventInfo.qScale();
616    
617     const gen::PdfInfo* pdf = genEventInfo.pdf();
618     if(pdf){
619       genInfo.pdf_id1=pdf->id.first;
620       genInfo.pdf_id2=pdf->id.second;
621       genInfo.pdf_x1=pdf->x.first;
622       genInfo.pdf_x2=pdf->x.second;
623       genInfo.pdf_scalePDF=pdf->scalePDF;
624       genInfo.pdf_xPDF1=pdf->xPDF.first;
625       genInfo.pdf_xPDF2=pdf->xPDF.second;
626     }
627     else{
628       genInfo.pdf_id1=-999;
629       genInfo.pdf_id2=-999;
630       genInfo.pdf_x1=-999;
631       genInfo.pdf_x2=-999;
632       genInfo.pdf_scalePDF=-999;
633       genInfo.pdf_xPDF1=-999;
634       genInfo.pdf_xPDF2=-999;
635     }
636
637     edm::Handle<std::vector<PileupSummaryInfo> > pus;
638     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
639     genInfo.pileup_NumInteractions_intime=0;
640     genInfo.pileup_NumInteractions_ootbefore=0;
641     genInfo.pileup_NumInteractions_ootafter=0;
642     if(pus.isValid()){
643       genInfo.pileup_TrueNumInteractions = (float) pus->at(0).getTrueNumInteractions();
644       for(size_t i=0; i<pus->size(); ++i){
645         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
646            genInfo.pileup_NumInteractions_intime += pus->at(i).getPU_NumInteractions();
647         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
648            genInfo.pileup_NumInteractions_ootbefore += pus->at(i).getPU_NumInteractions();
649         }
650         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
651           genInfo.pileup_NumInteractions_ootafter += pus->at(i).getPU_NumInteractions();
652         }
653       }
654     }
655
656     edm::Handle<reco::GenParticleCollection> genPartColl;
657     iEvent.getByLabel(edm::InputTag("genParticles"), genPartColl);
658     int index=-1;
659     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
660       index++;
661      
662       //write out only top quarks and status 3 particles (works fine only for MadGraph)
663       if(abs(iter->pdgId())==6 || iter->status()==3){
664         GenParticle genp;
665         genp.charge = iter->charge();;
666         genp.pt = iter->p4().pt();
667         genp.eta = iter->p4().eta();
668         genp.phi = iter->p4().phi();
669         genp.energy = iter->p4().E();
670         genp.index =index;
671         genp.status = iter->status();
672         genp.pdgId = iter->pdgId();
673
674         genp.mother1=-1;
675         genp.mother2=-1;
676         genp.daughter1=-1;
677         genp.daughter2=-1;
678        
679         int nm=iter->numberOfMothers();
680         int nd=iter->numberOfDaughters();
681        
682         if (nm>0) genp.mother1 = iter->motherRef(0).key();
683         if (nm>1) genp.mother2 = iter->motherRef(nm-1).key();
684         if (nd>0) genp.daughter1 = iter->daughterRef(0).key();
685         if (nd>1) genp.daughter2 = iter->daughterRef(nd-1).key();
686        
687         genps.push_back(genp);
688       }
689     }
690
691   }
692
857  
858     tr->Fill();
859 <
859 >   if(doLumiInfo)
860 >     previouslumiblockwasfilled=true;
861   }
862  
863  
# Line 700 | Line 865 | NtupleWriter::analyze(const edm::Event&
865   void
866   NtupleWriter::beginJob()
867   {
868 +  if(doLumiInfo){
869 +    totalRecLumi=0;
870 +    totalDelLumi=0;
871 +    previouslumiblockwasfilled=false;
872 +  }
873   }
874  
875   // ------------ method called once each job just after ending the event loop  ------------
# Line 716 | Line 886 | void
886   NtupleWriter::beginRun(edm::Run const& iRun, edm::EventSetup const&  iSetup)
887   {
888    if(doTrigger){
889 <    bool setup_changed = false;
890 <    hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
889 >    //bool setup_changed = false;
890 >    //hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
891      newrun=true;
892    }
893  
894    if(doJets || doTopJets){
895 <    edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
896 <    iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
897 <    JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
898 <    jecUnc = new JetCorrectionUncertainty(JetCorPar);
895 >    if(doJECUncertainty){
896 >      edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
897 >      iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
898 >      JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
899 >      jecUnc = new JetCorrectionUncertainty(JetCorPar);
900 >    }
901    }
902   }
903  
# Line 733 | Line 905 | NtupleWriter::beginRun(edm::Run const& i
905   void
906   NtupleWriter::endRun(edm::Run const&, edm::EventSetup const&)
907   {
908 +  if(doLumiInfo)
909 +    std::cout << "total integ. luminosity: " << totalDelLumi <<"(del) " << totalRecLumi << "(rec)" << std::endl;
910   }
911  
912   // ------------ method called when starting to processes a luminosity block  ------------
913   void
914 < NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&)
914 > NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const& lumi, edm::EventSetup const&)
915   {
916 +  if(doLumiInfo){
917 +    edm::Handle<LumiSummary> l;
918 +    lumi.getByLabel("lumiProducer", l);
919 +    
920 +    //add lumi of lumi blocks without any event to next lumiblock
921 +    if(previouslumiblockwasfilled){
922 +      intgRecLumi=0;
923 +      intgDelLumi=0;
924 +    }
925 +    previouslumiblockwasfilled=false;
926 +    
927 +    if (l.isValid()){;
928 +      intgRecLumi+=l->intgRecLumi()*6.37;
929 +      intgDelLumi+=l->intgDelLumi()*6.37;
930 +      totalRecLumi+=l->intgRecLumi()*6.37;
931 +      totalDelLumi+=l->intgDelLumi()*6.37;
932 +    }
933 +    //std::cout << "this lb: " <<l->intgRecLumi()*6.37 <<"   " << l->intgDelLumi()*6.37<<std::endl;
934 +    //std::cout << "summed:  "<< intgRecLumi << "   " << intgDelLumi << std::endl;
935 +  }
936   }
937  
938   // ------------ method called when ending the processing of a luminosity block  ------------

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines