ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc

Comparing UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc (file contents):
Revision 1.16 by peiffer, Mon Apr 30 09:01:16 2012 UTC vs.
Revision 1.26 by peiffer, Mon Jun 10 13:33:05 2013 UTC

#	Line 17 | Line 17
17		//
18		//
19
20	<	#include "UHHAnalysis//NtupleWriter/interface/NtupleWriter.h"
21	<
22	<
20	>	#include "UHHAnalysis/NtupleWriter/interface/NtupleWriter.h"
21	>	#include "RecoBTau/JetTagComputer/interface/GenericMVAJetTagComputer.h"
22	>	#include "RecoBTau/JetTagComputer/interface/GenericMVAJetTagComputerWrapper.h"
23	>	#include "RecoBTau/JetTagComputer/interface/JetTagComputer.h"
24	>	#include "RecoBTau/JetTagComputer/interface/JetTagComputerRecord.h"
25	>	#include "RecoBTag/SecondaryVertex/interface/CombinedSVComputer.h"
26
27		//
28		// constants, enums and typedefs
#	Line 50 | Line 53 | NtupleWriter::NtupleWriter(const edm::Pa
53		doMuons = iConfig.getParameter<bool>("doMuons");
54		doTaus = iConfig.getParameter<bool>("doTaus");
55		doJets = iConfig.getParameter<bool>("doJets");
56	+	doGenJets = iConfig.getParameter<bool>("doGenJets");
57	+	doJECUncertainty = iConfig.getParameter<bool>("doJECUncertainty");
58		doGenTopJets = iConfig.getParameter<bool>("doGenTopJets");
59		doPhotons = iConfig.getParameter<bool>("doPhotons");
60		doMET = iConfig.getParameter<bool>("doMET");
#	Line 59 | Line 64 | NtupleWriter::NtupleWriter(const edm::Pa
64		doPV = iConfig.getParameter<bool>("doPV");
65		doTopJets = iConfig.getParameter<bool>("doTopJets");
66		doTrigger = iConfig.getParameter<bool>("doTrigger");
67	<
67	>	SVComputer_  = iConfig.getUntrackedParameter<edm::InputTag>("svComputer",edm::InputTag("combinedSecondaryVertex"));
68	>	doTagInfos = iConfig.getUntrackedParameter<bool>("doTagInfos",false);
69		// initialization of tree variables
70
71		tr->Branch("run",&run);
72		tr->Branch("event",&event);
73		tr->Branch("luminosityBlock",&luminosityBlock);
74		tr->Branch("isRealData",&isRealData);
75	<	tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
75	>	tr->Branch("rho",&rho);
76	>	rho_source = iConfig.getParameter<edm::InputTag>("rho_source");
77	>
78	>	//tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
79		if(doLumiInfo){
80		tr->Branch("intgRecLumi",&intgRecLumi);
81		tr->Branch("intgDelLumi",&intgDelLumi);
#	Line 104 | Line 113 | NtupleWriter::NtupleWriter(const edm::Pa
113		tr->Branch( jet_sources[j].c_str(), "std::vector<Jet>", &jets[j]);
114		}
115		}
116	+	if(doGenJets){
117	+	genjet_sources = iConfig.getParameter<std::vector<std::string> >("genjet_sources");
118	+	genjet_ptmin = iConfig.getParameter<double> ("genjet_ptmin");
119	+	genjet_etamax = iConfig.getParameter<double> ("genjet_etamax");
120	+	for(size_t j=0; j< genjet_sources.size(); ++j){
121	+	tr->Branch( genjet_sources[j].c_str(), "std::vector<Particle>", &genjets[j]);
122	+	}
123	+	}
124		if(doTopJets){
125		topjet_sources = iConfig.getParameter<std::vector<std::string> >("topjet_sources");
126		topjet_ptmin = iConfig.getParameter<double> ("topjet_ptmin");
#	Line 139 | Line 156 | NtupleWriter::NtupleWriter(const edm::Pa
156		}
157		}
158		if(doGenInfo){
159	+	genparticle_source= iConfig.getParameter<edm::InputTag>("genparticle_source");
160		tr->Branch("genInfo","GenInfo",&genInfo);
161		tr->Branch("GenParticles","std::vector<GenParticle>", &genps);
162		}
#	Line 181 | Line 199 | NtupleWriter::analyze(const edm::Event&
199		luminosityBlock = iEvent.luminosityBlock();
200		isRealData      = iEvent.isRealData();
201
202	<	if(isRealData){
203	<	edm::Handle<bool> bool_handle;
204	<	iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
205	<	HBHENoiseFilterResult = *(bool_handle.product());
206	<	}
207	<	else HBHENoiseFilterResult = false;
202	>	edm::Handle<double> m_rho;
203	>	iEvent.getByLabel(rho_source,m_rho);
204	>	rho=*m_rho;
205	>
206	>	//    if(isRealData){
207	>	//      edm::Handle<bool> bool_handle;
208	>	//      iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
209	>	//      HBHENoiseFilterResult = *(bool_handle.product());
210	>	//    }
211	>	//    else HBHENoiseFilterResult = false;
212
213		// ------------- primary vertices and beamspot  -------------
214
#	Line 202 | Line 224 | NtupleWriter::analyze(const edm::Event&
224		reco::Vertex reco_pv = reco_pvs[i];
225
226		PrimaryVertex pv;
227	<	pv.x =  reco_pv.x();
228	<	pv.y =  reco_pv.y();
229	<	pv.z =  reco_pv.z();
230	<	pv.nTracks =  reco_pv.nTracks();
231	<	//pv.isValid =  reco_pv.isValid();
232	<	pv.chi2 =  reco_pv.chi2();
233	<	pv.ndof =  reco_pv.ndof();
227	>	pv.set_x( reco_pv.x());
228	>	pv.set_y( reco_pv.y());
229	>	pv.set_z( reco_pv.z());
230	>	pv.set_nTracks( reco_pv.nTracks());
231	>	//pv.set_isValid( reco_pv.isValid());
232	>	pv.set_chi2( reco_pv.chi2());
233	>	pv.set_ndof( reco_pv.ndof());
234
235		pvs[j].push_back(pv);
236		}
#	Line 226 | Line 248 | NtupleWriter::analyze(const edm::Event&
248		// ------------- generator info -------------
249
250		if(doGenInfo){
251	<	genInfo.weights.clear();
252	<	genInfo.binningValues.clear();
251	>	genInfo.clear_weights();
252	>	genInfo.clear_binningValues();
253		genps.clear();
254
255		edm::Handle<GenEventInfoProduct> genEventInfoProduct;
256		iEvent.getByLabel("generator", genEventInfoProduct);
257		const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
258
259	<	genInfo.binningValues = genEventInfo.binningValues();
260	<	genInfo.weights = genEventInfo.weights();
261	<	genInfo.alphaQCD = genEventInfo.alphaQCD();
262	<	genInfo.alphaQED = genEventInfo.alphaQED();
263	<	genInfo.qScale = genEventInfo.qScale();
259	>	for(unsigned int k=0; k<genEventInfo.binningValues().size();++k){
260	>	genInfo.add_binningValue(genEventInfo.binningValues().at(k));
261	>	}
262	>	for(unsigned int k=0; k<genEventInfo.weights().size();++k){
263	>	genInfo.add_weight(genEventInfo.weights().at(k));
264	>	}
265	>	genInfo.set_alphaQCD(genEventInfo.alphaQCD());
266	>	genInfo.set_alphaQED(genEventInfo.alphaQED());
267	>	genInfo.set_qScale(genEventInfo.qScale());
268
269		const gen::PdfInfo* pdf = genEventInfo.pdf();
270		if(pdf){
271	<	genInfo.pdf_id1=pdf->id.first;
272	<	genInfo.pdf_id2=pdf->id.second;
273	<	genInfo.pdf_x1=pdf->x.first;
274	<	genInfo.pdf_x2=pdf->x.second;
275	<	genInfo.pdf_scalePDF=pdf->scalePDF;
276	<	genInfo.pdf_xPDF1=pdf->xPDF.first;
277	<	genInfo.pdf_xPDF2=pdf->xPDF.second;
271	>	genInfo.set_pdf_id1(pdf->id.first);
272	>	genInfo.set_pdf_id2(pdf->id.second);
273	>	genInfo.set_pdf_x1(pdf->x.first);
274	>	genInfo.set_pdf_x2(pdf->x.second);
275	>	genInfo.set_pdf_scalePDF(pdf->scalePDF);
276	>	genInfo.set_pdf_xPDF1(pdf->xPDF.first);
277	>	genInfo.set_pdf_xPDF2(pdf->xPDF.second);
278		}
279		else{
280	<	genInfo.pdf_id1=-999;
281	<	genInfo.pdf_id2=-999;
282	<	genInfo.pdf_x1=-999;
283	<	genInfo.pdf_x2=-999;
284	<	genInfo.pdf_scalePDF=-999;
285	<	genInfo.pdf_xPDF1=-999;
286	<	genInfo.pdf_xPDF2=-999;
280	>	genInfo.set_pdf_id1(-999);
281	>	genInfo.set_pdf_id2(-999);
282	>	genInfo.set_pdf_x1(-999);
283	>	genInfo.set_pdf_x2(-999);
284	>	genInfo.set_pdf_scalePDF(-999);
285	>	genInfo.set_pdf_xPDF1(-999);
286	>	genInfo.set_pdf_xPDF2(-999);
287		}
288
289		edm::Handle<std::vector<PileupSummaryInfo> > pus;
290		iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
291	<	genInfo.pileup_NumInteractions_intime=0;
292	<	genInfo.pileup_NumInteractions_ootbefore=0;
293	<	genInfo.pileup_NumInteractions_ootafter=0;
291	>	genInfo.set_pileup_NumInteractions_intime(0);
292	>	genInfo.set_pileup_NumInteractions_ootbefore(0);
293	>	genInfo.set_pileup_NumInteractions_ootafter(0);
294		if(pus.isValid()){
295	<	genInfo.pileup_TrueNumInteractions = (float) pus->at(0).getTrueNumInteractions();
295	>	genInfo.set_pileup_TrueNumInteractions ( (float) pus->at(0).getTrueNumInteractions());
296		for(size_t i=0; i<pus->size(); ++i){
297		if(pus->at(i).getBunchCrossing() == 0) // intime pileup
298	<	genInfo.pileup_NumInteractions_intime += pus->at(i).getPU_NumInteractions();
298	>	genInfo.set_pileup_NumInteractions_intime( genInfo.pileup_NumInteractions_intime() + pus->at(i).getPU_NumInteractions());
299		else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
300	<	genInfo.pileup_NumInteractions_ootbefore += pus->at(i).getPU_NumInteractions();
300	>	genInfo.set_pileup_NumInteractions_ootbefore( genInfo.pileup_NumInteractions_ootbefore() + pus->at(i).getPU_NumInteractions());
301		}
302		else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
303	<	genInfo.pileup_NumInteractions_ootafter += pus->at(i).getPU_NumInteractions();
303	>	genInfo.set_pileup_NumInteractions_ootafter( genInfo.pileup_NumInteractions_ootafter() + pus->at(i).getPU_NumInteractions());
304		}
305		}
306		}
307
308		edm::Handle<reco::GenParticleCollection> genPartColl;
309	<	iEvent.getByLabel(edm::InputTag("genParticles"), genPartColl);
309	>	iEvent.getByLabel(genparticle_source, genPartColl);
310		int index=-1;
311		for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
312		index++;
313
314	<	//write out only top quarks and status 3 particles (works fine only for MadGraph)
315	<	if(abs(iter->pdgId())==6 || iter->status()==3 || doAllGenParticles){
314	>	//write out only top quarks,final state leptons and status 3 particles (works fine only for MadGraph)
315	>	bool islepton = iter->status()==1 && abs(iter->pdgId())>=11 && abs(iter->pdgId())<=16 ;
316	>	if(abs(iter->pdgId())==6 || iter->status()==3 || islepton ||  doAllGenParticles){
317		GenParticle genp;
318	<	genp.charge = iter->charge();;
319	<	genp.pt = iter->p4().pt();
320	<	genp.eta = iter->p4().eta();
321	<	genp.phi = iter->p4().phi();
322	<	genp.energy = iter->p4().E();
323	<	genp.index =index;
324	<	genp.status = iter->status();
325	<	genp.pdgId = iter->pdgId();
326	<
327	<	genp.mother1=-1;
328	<	genp.mother2=-1;
329	<	genp.daughter1=-1;
330	<	genp.daughter2=-1;
318	>	genp.set_charge(iter->charge());
319	>	genp.set_pt(iter->p4().pt());
320	>	genp.set_eta(iter->p4().eta());
321	>	genp.set_phi(iter->p4().phi());
322	>	genp.set_energy(iter->p4().E());
323	>	genp.set_index(index);
324	>	genp.set_status( iter->status());
325	>	genp.set_pdgId( iter->pdgId());
326	>
327	>	genp.set_mother1(-1);
328	>	genp.set_mother2(-1);
329	>	genp.set_daughter1(-1);
330	>	genp.set_daughter2(-1);
331
332		int nm=iter->numberOfMothers();
333		int nd=iter->numberOfDaughters();
334
335
336	<	if (nm>0) genp.mother1 = iter->motherRef(0).key();
337	<	if (nm>1) genp.mother2 = iter->motherRef(1).key();
338	<	if (nd>0) genp.daughter1 = iter->daughterRef(0).key();
339	<	if (nd>1) genp.daughter2 = iter->daughterRef(1).key();
340	<	//std::cout << genp.index <<"  pdgId = " << genp.pdgId << "  mo1 = " << genp.mother1 << "  mo2 = " << genp.mother2 <<"  da1 = " << genp.daughter1 << "  da2 = " << genp.daughter2 <<std::endl;
336	>	if (nm>0) genp.set_mother1( iter->motherRef(0).key());
337	>	if (nm>1) genp.set_mother2( iter->motherRef(1).key());
338	>	if (nd>0) genp.set_daughter1( iter->daughterRef(0).key());
339	>	if (nd>1) genp.set_daughter2( iter->daughterRef(1).key());
340	>
341		genps.push_back(genp);
342		}
343		}
#	Line 319 | Line 346 | NtupleWriter::analyze(const edm::Event&
346		// ------------- electrons -------------
347		if(doElectrons){
348
349	<	edm::Handle<reco::ConversionCollection> hConversions;
350	<	iEvent.getByLabel("allConversions", hConversions);
349	>	//      edm::Handle<reco::ConversionCollection> hConversions;
350	>	//      iEvent.getByLabel("allConversions", hConversions);
351
352	<	edm::Handle<reco::BeamSpot> beamSpot;
353	<	iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
354	<	const reco::BeamSpot & bsp = *beamSpot;
352	>	//      edm::Handle<reco::BeamSpot> beamSpot;
353	>	//      iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
354	>	//      const reco::BeamSpot & bsp = *beamSpot;
355
356		for(size_t j=0; j< electron_sources.size(); ++j){
357		eles[j].clear();
#	Line 336 | Line 363 | NtupleWriter::analyze(const edm::Event&
363		pat::Electron pat_ele = pat_electrons[i];
364		Electron ele;
365
366	<	ele.charge =  pat_ele.charge();
367	<	ele.pt =  pat_ele.pt();
368	<	ele.eta =  pat_ele.eta();
369	<	ele.phi =  pat_ele.phi();
370	<	ele.energy =  pat_ele.energy();
371	<	ele.vertex_x = pat_ele.vertex().x();
372	<	ele.vertex_y = pat_ele.vertex().y();
373	<	ele.vertex_z = pat_ele.vertex().z();
374	<	ele.supercluster_eta = pat_ele.superCluster()->eta();
375	<	ele.supercluster_phi = pat_ele.superCluster()->phi();
376	<	ele.dB = pat_ele.dB();
377	<	//ele.particleIso = pat_ele.particleIso();
378	<	ele.neutralHadronIso = pat_ele.neutralHadronIso();
379	<	ele.chargedHadronIso = pat_ele.chargedHadronIso();
380	<	ele.trackIso = pat_ele.trackIso();
381	<	ele.photonIso = pat_ele.photonIso();
382	<	ele.puChargedHadronIso = pat_ele.puChargedHadronIso();
383	<	ele.gsfTrack_trackerExpectedHitsInner_numberOfLostHits = pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits();
384	<	ele.gsfTrack_px= pat_ele.gsfTrack()->px();
385	<	ele.gsfTrack_py= pat_ele.gsfTrack()->py();
386	<	ele.gsfTrack_pz= pat_ele.gsfTrack()->pz();
387	<	ele.gsfTrack_vx= pat_ele.gsfTrack()->vx();
388	<	ele.gsfTrack_vy= pat_ele.gsfTrack()->vy();
389	<	ele.gsfTrack_vz= pat_ele.gsfTrack()->vz();
390	<	ele.passconversionveto = !ConversionTools::hasMatchedConversion(pat_ele,hConversions,bsp.position());
391	<	ele.dEtaIn=pat_ele.deltaEtaSuperClusterTrackAtVtx();
392	<	ele.dPhiIn=pat_ele.deltaPhiSuperClusterTrackAtVtx();
393	<	ele.sigmaIEtaIEta=pat_ele.sigmaIetaIeta();
394	<	ele.HoverE=pat_ele.hadronicOverEm();
395	<	ele.fbrem=pat_ele.fbrem();
396	<	ele.EoverPIn=pat_ele.eSuperClusterOverP();
397	<	ele.EcalEnergy=pat_ele.ecalEnergy();
398	<	//ele.mvaTrigV0=pat_ele.electronID("mvaTrigV0");
399	<	//ele.mvaNonTrigV0=pat_ele.electronID("mvaNonTrigV0");
366	>	ele.set_charge( pat_ele.charge());
367	>	ele.set_pt( pat_ele.pt());
368	>	ele.set_eta( pat_ele.eta());
369	>	ele.set_phi( pat_ele.phi());
370	>	ele.set_energy( pat_ele.energy());
371	>	ele.set_vertex_x(pat_ele.vertex().x());
372	>	ele.set_vertex_y(pat_ele.vertex().y());
373	>	ele.set_vertex_z(pat_ele.vertex().z());
374	>	ele.set_supercluster_eta(pat_ele.superCluster()->eta());
375	>	ele.set_supercluster_phi(pat_ele.superCluster()->phi());
376	>	ele.set_dB(pat_ele.dB());
377	>	//ele.set_particleIso(pat_ele.particleIso());
378	>	ele.set_neutralHadronIso(pat_ele.neutralHadronIso());
379	>	ele.set_chargedHadronIso(pat_ele.chargedHadronIso());
380	>	ele.set_trackIso(pat_ele.trackIso());
381	>	ele.set_photonIso(pat_ele.photonIso());
382	>	ele.set_puChargedHadronIso(pat_ele.puChargedHadronIso());
383	>	ele.set_gsfTrack_trackerExpectedHitsInner_numberOfLostHits(pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits());
384	>	ele.set_gsfTrack_px( pat_ele.gsfTrack()->px());
385	>	ele.set_gsfTrack_py( pat_ele.gsfTrack()->py());
386	>	ele.set_gsfTrack_pz( pat_ele.gsfTrack()->pz());
387	>	ele.set_gsfTrack_vx( pat_ele.gsfTrack()->vx());
388	>	ele.set_gsfTrack_vy( pat_ele.gsfTrack()->vy());
389	>	ele.set_gsfTrack_vz( pat_ele.gsfTrack()->vz());
390	>	//ele.set_passconversionveto(!ConversionTools::hasMatchedConversion(pat_ele,hConversions,bsp.position()));
391	>	ele.set_passconversionveto(pat_ele.passConversionVeto());
392	>	ele.set_dEtaIn(pat_ele.deltaEtaSuperClusterTrackAtVtx());
393	>	ele.set_dPhiIn(pat_ele.deltaPhiSuperClusterTrackAtVtx());
394	>	ele.set_sigmaIEtaIEta(pat_ele.sigmaIetaIeta());
395	>	ele.set_HoverE(pat_ele.hadronicOverEm());
396	>	ele.set_fbrem(pat_ele.fbrem());
397	>	ele.set_EoverPIn(pat_ele.eSuperClusterOverP());
398	>	ele.set_EcalEnergy(pat_ele.ecalEnergy());
399	>	ele.set_mvaTrigV0(pat_ele.electronID("mvaTrigV0"));
400	>	ele.set_mvaNonTrigV0(pat_ele.electronID("mvaNonTrigV0"));
401	>	float AEff03 = 0.00;
402	>	if(isRealData){
403	>	AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAData2011);
404	>	}else{
405	>	AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAFall11MC);
406	>	}
407	>	ele.set_AEff(AEff03);
408
409		eles[j].push_back(ele);
410		}
#	Line 389 | Line 424 | NtupleWriter::analyze(const edm::Event&
424		pat::Muon pat_mu = pat_muons[i];
425
426		Muon mu;
427	<	mu.charge =  pat_mu.charge();
428	<	mu.pt =  pat_mu.pt();
429	<	mu.eta =  pat_mu.eta();
430	<	mu.phi =  pat_mu.phi();
431	<	mu.energy =  pat_mu.energy();
432	<	mu.vertex_x = pat_mu.vertex().x();
433	<	mu.vertex_y = pat_mu.vertex().y();
434	<	mu.vertex_z = pat_mu.vertex().z();
435	<	mu.dB = pat_mu.dB();
436	<	//mu.particleIso = pat_mu.particleIso();
437	<	mu.neutralHadronIso = pat_mu.neutralHadronIso();
438	<	mu.chargedHadronIso = pat_mu.chargedHadronIso();
439	<	mu.trackIso = pat_mu.trackIso();
440	<	mu.photonIso = pat_mu.photonIso();
441	<	mu.puChargedHadronIso = pat_mu.puChargedHadronIso();
442	<	mu.isGlobalMuon = pat_mu.isGlobalMuon();
443	<	mu.isStandAloneMuon = pat_mu.isStandAloneMuon();
444	<	mu.isTrackerMuon = pat_mu.isTrackerMuon();
445	<	mu.numberOfMatchedStations = pat_mu.numberOfMatchedStations();
427	>	mu.set_charge( pat_mu.charge());
428	>	mu.set_pt( pat_mu.pt());
429	>	mu.set_eta( pat_mu.eta());
430	>	mu.set_phi( pat_mu.phi());
431	>	mu.set_energy( pat_mu.energy());
432	>	mu.set_vertex_x ( pat_mu.vertex().x());
433	>	mu.set_vertex_y ( pat_mu.vertex().y());
434	>	mu.set_vertex_z ( pat_mu.vertex().z());
435	>	mu.set_dB ( pat_mu.dB());
436	>	//mu.particleIso ( pat_mu.particleIso());
437	>	mu.set_neutralHadronIso ( pat_mu.neutralHadronIso());
438	>	mu.set_chargedHadronIso ( pat_mu.chargedHadronIso());
439	>	mu.set_trackIso ( pat_mu.trackIso());
440	>	mu.set_photonIso ( pat_mu.photonIso());
441	>	mu.set_puChargedHadronIso ( pat_mu.puChargedHadronIso());
442	>	mu.set_isGlobalMuon ( pat_mu.isGlobalMuon());
443	>	mu.set_isPFMuon ( pat_mu.isPFMuon());
444	>	mu.set_isStandAloneMuon ( pat_mu.isStandAloneMuon());
445	>	mu.set_isTrackerMuon ( pat_mu.isTrackerMuon());
446	>	mu.set_numberOfMatchedStations ( pat_mu.numberOfMatchedStations());
447		reco::TrackRef globalTrack = pat_mu.globalTrack();
448		if(!globalTrack.isNull()){
449	<	mu.globalTrack_chi2 = globalTrack->chi2();
450	<	mu.globalTrack_ndof = globalTrack->ndof();
451	<	mu.globalTrack_d0 = globalTrack->d0();
452	<	mu.globalTrack_d0Error = globalTrack->d0Error();
453	<	mu.globalTrack_numberOfValidHits = globalTrack->numberOfValidHits();
454	<	mu.globalTrack_numberOfLostHits = globalTrack->numberOfLostHits();
449	>	mu.set_globalTrack_chi2 ( globalTrack->chi2());
450	>	mu.set_globalTrack_ndof ( globalTrack->ndof());
451	>	mu.set_globalTrack_d0 ( globalTrack->d0());
452	>	mu.set_globalTrack_d0Error ( globalTrack->d0Error());
453	>	mu.set_globalTrack_numberOfValidHits ( globalTrack->numberOfValidHits());
454	>	mu.set_globalTrack_numberOfLostHits ( globalTrack->numberOfLostHits());
455	>	mu.set_globalTrack_numberOfValidMuonHits(globalTrack->hitPattern().numberOfValidMuonHits() );
456		}
457		else{
458	<	mu.globalTrack_chi2 = 0;
459	<	mu.globalTrack_ndof = 0;
460	<	mu.globalTrack_d0 = 0;
461	<	mu.globalTrack_d0Error = 0;
462	<	mu.globalTrack_numberOfValidHits = 0;
463	<	mu.globalTrack_numberOfLostHits = 0;
458	>	mu.set_globalTrack_chi2 ( 0);
459	>	mu.set_globalTrack_ndof ( 0);
460	>	mu.set_globalTrack_d0 ( 0);
461	>	mu.set_globalTrack_d0Error ( 0);
462	>	mu.set_globalTrack_numberOfValidHits ( 0);
463	>	mu.set_globalTrack_numberOfLostHits ( 0);
464		}
465		reco::TrackRef innerTrack = pat_mu.innerTrack();
466		if(!innerTrack.isNull()){
467	<	mu.innerTrack_chi2 = innerTrack->chi2();
468	<	mu.innerTrack_ndof = innerTrack->ndof();
469	<	mu.innerTrack_d0 = innerTrack->d0();
470	<	mu.innerTrack_d0Error = innerTrack->d0Error();
471	<	mu.innerTrack_numberOfValidHits = innerTrack->numberOfValidHits();
472	<	mu.innerTrack_numberOfLostHits = innerTrack->numberOfLostHits();
473	<	mu.innerTrack_trackerLayersWithMeasurement = innerTrack->hitPattern().trackerLayersWithMeasurement();
474	<	mu.innerTrack_numberOfValidPixelHits = innerTrack->hitPattern().numberOfValidPixelHits();
467	>	mu.set_innerTrack_chi2 ( innerTrack->chi2());
468	>	mu.set_innerTrack_ndof ( innerTrack->ndof());
469	>	mu.set_innerTrack_d0 ( innerTrack->d0());
470	>	mu.set_innerTrack_d0Error ( innerTrack->d0Error());
471	>	mu.set_innerTrack_numberOfValidHits ( innerTrack->numberOfValidHits());
472	>	mu.set_innerTrack_numberOfLostHits ( innerTrack->numberOfLostHits());
473	>	mu.set_innerTrack_trackerLayersWithMeasurement ( innerTrack->hitPattern().trackerLayersWithMeasurement());
474	>	mu.set_innerTrack_numberOfValidPixelHits ( innerTrack->hitPattern().numberOfValidPixelHits());
475		}
476		else{
477	<	mu.innerTrack_chi2 = 0;
478	<	mu.innerTrack_ndof = 0;
479	<	mu.innerTrack_d0 = 0;
480	<	mu.innerTrack_d0Error = 0;
481	<	mu.innerTrack_numberOfValidHits = 0;
482	<	mu.innerTrack_numberOfLostHits = 0;
483	<	mu.innerTrack_trackerLayersWithMeasurement = 0;
484	<	mu.innerTrack_numberOfValidPixelHits = 0;
477	>	mu.set_innerTrack_chi2 ( 0);
478	>	mu.set_innerTrack_ndof ( 0);
479	>	mu.set_innerTrack_d0 ( 0);
480	>	mu.set_innerTrack_d0Error ( 0);
481	>	mu.set_innerTrack_numberOfValidHits ( 0);
482	>	mu.set_innerTrack_numberOfLostHits ( 0);
483	>	mu.set_innerTrack_trackerLayersWithMeasurement ( 0);
484	>	mu.set_innerTrack_numberOfValidPixelHits ( 0);
485		}
486		reco::TrackRef outerTrack = pat_mu.outerTrack();
487		if(!outerTrack.isNull()){
488	<	mu.outerTrack_chi2 = outerTrack->chi2();
489	<	mu.outerTrack_ndof = outerTrack->ndof();
490	<	mu.outerTrack_d0 = outerTrack->d0();
491	<	mu.outerTrack_d0Error = outerTrack->d0Error();
492	<	mu.outerTrack_numberOfValidHits = outerTrack->numberOfValidHits();
493	<	mu.outerTrack_numberOfLostHits = outerTrack->numberOfLostHits();
488	>	mu.set_outerTrack_chi2 ( outerTrack->chi2());
489	>	mu.set_outerTrack_ndof ( outerTrack->ndof());
490	>	mu.set_outerTrack_d0 ( outerTrack->d0());
491	>	mu.set_outerTrack_d0Error ( outerTrack->d0Error());
492	>	mu.set_outerTrack_numberOfValidHits ( outerTrack->numberOfValidHits());
493	>	mu.set_outerTrack_numberOfLostHits ( outerTrack->numberOfLostHits());
494		}
495		else{
496	<	mu.outerTrack_chi2 = 0;
497	<	mu.outerTrack_ndof = 0;
498	<	mu.outerTrack_d0 = 0;
499	<	mu.outerTrack_d0Error = 0;
500	<	mu.outerTrack_numberOfValidHits = 0;
501	<	mu.outerTrack_numberOfLostHits = 0;
496	>	mu.set_outerTrack_chi2 ( 0);
497	>	mu.set_outerTrack_ndof ( 0);
498	>	mu.set_outerTrack_d0 ( 0);
499	>	mu.set_outerTrack_d0Error ( 0);
500	>	mu.set_outerTrack_numberOfValidHits ( 0);
501	>	mu.set_outerTrack_numberOfLostHits ( 0);
502		}
503
504		mus[j].push_back(mu);
#	Line 485 | Line 522 | NtupleWriter::analyze(const edm::Event&
522		if(fabs(pat_tau.eta()) > tau_etamax) continue;
523
524		Tau tau;
525	<	tau.charge =  pat_tau.charge();
526	<	tau.pt =  pat_tau.pt();
527	<	tau.eta =  pat_tau.eta();
528	<	tau.phi =  pat_tau.phi();
529	<	tau.energy =  pat_tau.energy();
530	<	tau.decayModeFinding = pat_tau.tauID("decayModeFinding")>0.5;
531	<	tau.byVLooseCombinedIsolationDeltaBetaCorr  = pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5;
532	<	tau.byLooseCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5;
533	<	tau.byMediumCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5;
534	<	tau.byTightCombinedIsolationDeltaBetaCorr = pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5;
535	<	tau.againstElectronLoose  = pat_tau.tauID("againstElectronLoose")>0.5;
536	<	tau.againstElectronMedium = pat_tau.tauID("againstElectronMedium")>0.5;
537	<	tau.againstElectronTight = pat_tau.tauID("againstElectronTight")>0.5;
538	<	tau.againstElectronMVA  = pat_tau.tauID("againstElectronMVA")>0.5;
539	<	tau.againstMuonLoose = pat_tau.tauID("againstMuonLoose")>0.5;
540	<	tau.againstMuonMedium = pat_tau.tauID("againstMuonMedium")>0.5;
541	<	tau.againstMuonTight = pat_tau.tauID("againstMuonTight")>0.5;
542	<
543	<	reco::PFCandidateRef leadPFCand = pat_tau.leadPFCand();
544	<	if(!leadPFCand.isNull()){
545	<	tau.leadPFCand_px = leadPFCand->px();
546	<	tau.leadPFCand_py = leadPFCand->py();
547	<	tau.leadPFCand_pz = leadPFCand->pz();
548	<	}
549	<	else{
550	<	tau.leadPFCand_px = 0;
551	<	tau.leadPFCand_py = 0;
552	<	tau.leadPFCand_pz = 0;
553	<	}
525	>	tau.set_charge( pat_tau.charge());
526	>	tau.set_pt( pat_tau.pt());
527	>	tau.set_eta( pat_tau.eta());
528	>	tau.set_phi( pat_tau.phi());
529	>	tau.set_energy( pat_tau.energy());
530	>	tau.set_decayModeFinding ( pat_tau.tauID("decayModeFinding")>0.5);
531	>	//tau.set_byVLooseCombinedIsolationDeltaBetaCorr  ( pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5);
532	>	tau.set_byLooseCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5);
533	>	tau.set_byMediumCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5);
534	>	tau.set_byTightCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5);
535	>	tau.set_byLooseIsolationMVA( pat_tau.tauID("byLooseIsolationMVA")>0.5);
536	>	tau.set_byMediumIsolationMVA( pat_tau.tauID("byMediumIsolationMVA")>0.5);
537	>	tau.set_byTightIsolationMVA( pat_tau.tauID("byTightIsolationMVA")>0.5);
538	>	tau.set_byLooseIsolationMVA2( pat_tau.tauID("byLooseIsolationMVA2")>0.5);
539	>	tau.set_byMediumIsolationMVA2( pat_tau.tauID("byMediumIsolationMVA2")>0.5);
540	>	tau.set_byTightIsolationMVA2( pat_tau.tauID("byTightIsolationMVA2")>0.5);
541	>	tau.set_byLooseCombinedIsolationDeltaBetaCorr3Hits(  pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr3Hits")>0.5);
542	>	tau.set_byMediumCombinedIsolationDeltaBetaCorr3Hits ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr3Hits")>0.5);
543	>	tau.set_byTightCombinedIsolationDeltaBetaCorr3Hits ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr3Hits")>0.5);
544	>	tau.set_againstElectronLooseMVA3  ( pat_tau.tauID("againstElectronLooseMVA3")>0.5);
545	>	tau.set_againstElectronMediumMVA3 ( pat_tau.tauID("againstElectronMediumMVA3")>0.5);
546	>	tau.set_againstElectronTightMVA3 ( pat_tau.tauID("againstElectronTightMVA3")>0.5);
547	>	tau.set_againstElectronVTightMVA3 ( pat_tau.tauID("againstElectronVTightMVA3")>0.5);
548	>	tau.set_againstMuonLoose2 ( pat_tau.tauID("againstMuonLoose2")>0.5);
549	>	tau.set_againstMuonMedium2 ( pat_tau.tauID("againstMuonMedium2")>0.5);
550	>	tau.set_againstMuonTight2 ( pat_tau.tauID("againstMuonTight2")>0.5);
551	>	tau.set_byIsolationMVAraw(  pat_tau.tauID("byIsolationMVAraw"));
552	>	tau.set_byIsolationMVA2raw(  pat_tau.tauID("byIsolationMVA2raw"));
553	>	tau.set_decayMode( pat_tau.decayMode() );
554	>	tau.set_byCombinedIsolationDeltaBetaCorrRaw( pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw"));
555	>	tau.set_byCombinedIsolationDeltaBetaCorrRaw3Hits( pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw3Hits"));
556	>
557	>	//       std::cout << pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw3Hits") << std::endl;
558	>
559	>	//       reco::PFCandidateRef leadPFCand = pat_tau.leadPFCand();
560	>	//       if(!leadPFCand.isNull()){
561	>	//         tau.set_leadPFCand_px ( leadPFCand->px());
562	>	//         tau.set_leadPFCand_py ( leadPFCand->py());
563	>	//         tau.set_leadPFCand_pz ( leadPFCand->pz());
564	>	//       }
565	>	//       else{
566	>	//         tau.set_leadPFCand_px ( 0);
567	>	//         tau.set_leadPFCand_py ( 0);
568	>	//         tau.set_leadPFCand_pz ( 0);
569	>	//       }
570		taus[j].push_back(tau);
571		}
572		}
573		}
574
575	+	//-------------- gen jets -------------
576	+
577	+	if(doGenJets){
578	+	for(size_t j=0; j< genjet_sources.size(); ++j){
579	+
580	+	genjets[j].clear();
581	+
582	+	edm::Handle< std::vector<reco::GenJet> > genjet_handle;
583	+	iEvent.getByLabel(genjet_sources[j], genjet_handle);
584	+	const std::vector<reco::GenJet>& gen_jets = *(genjet_handle.product());
585	+
586	+	for (unsigned int i = 0; i < gen_jets.size(); ++i) {
587	+	const reco::GenJet* gen_jet = &gen_jets[i];
588	+	if(gen_jet->pt() < genjet_ptmin) continue;
589	+	if(fabs(gen_jet->eta()) > genjet_etamax) continue;
590	+
591	+	Particle jet;
592	+	jet.set_charge(gen_jet->charge());
593	+	jet.set_pt(gen_jet->pt());
594	+	jet.set_eta(gen_jet->eta());
595	+	jet.set_phi(gen_jet->phi());
596	+	jet.set_energy(gen_jet->energy());
597	+
598	+	genjets[j].push_back(jet);
599	+
600	+	}
601	+	}
602	+	}
603	+
604		// ------------- jets -------------
605		if(doJets){
606		for(size_t j=0; j< jet_sources.size(); ++j){
#	Line 540 | Line 622 | NtupleWriter::analyze(const edm::Event&
622		//       }
623
624		Jet jet;
625	<	jet.charge = pat_jet.charge();
626	<	jet.pt = pat_jet.pt();
627	<	jet.eta = pat_jet.eta();
628	<	jet.phi = pat_jet.phi();
629	<	jet.energy = pat_jet.energy();
630	<	jet.numberOfDaughters =pat_jet.numberOfDaughters();
625	>	jet.set_charge(pat_jet.charge());
626	>	jet.set_pt(pat_jet.pt());
627	>	jet.set_eta(pat_jet.eta());
628	>	jet.set_phi(pat_jet.phi());
629	>	jet.set_energy(pat_jet.energy());
630	>	jet.set_numberOfDaughters (pat_jet.numberOfDaughters());
631		const reco::TrackRefVector&  jettracks = pat_jet.associatedTracks();
632	<	jet.nTracks = jettracks.size();
633	<	jet.jetArea = pat_jet.jetArea();
634	<	jet.pileup = pat_jet.pileup();
632	>	jet.set_nTracks ( jettracks.size());
633	>	jet.set_jetArea(pat_jet.jetArea());
634	>	jet.set_pileup(pat_jet.pileup());
635		if(pat_jet.isPFJet()){
636	<	jet.neutralEmEnergyFraction =pat_jet.neutralEmEnergyFraction();
637	<	jet.neutralHadronEnergyFraction =pat_jet.neutralHadronEnergyFraction();
638	<	jet.chargedEmEnergyFraction =pat_jet.chargedEmEnergyFraction();
639	<	jet.chargedHadronEnergyFraction =pat_jet.chargedHadronEnergyFraction();
640	<	jet.muonEnergyFraction =pat_jet.muonEnergyFraction();
641	<	jet.photonEnergyFraction =pat_jet.photonEnergyFraction();
642	<	jet.chargedMultiplicity =pat_jet.chargedMultiplicity();
643	<	jet.neutralMultiplicity =pat_jet.neutralMultiplicity();
644	<	jet.muonMultiplicity =pat_jet.muonMultiplicity();
645	<	jet.electronMultiplicity =pat_jet.electronMultiplicity();
646	<	jet.photonMultiplicity =pat_jet.photonMultiplicity();
647	<	}
648	<
649	<	jecUnc->setJetEta(pat_jet.eta());
650	<	jecUnc->setJetPt(pat_jet.pt());
651	<	jet.JEC_uncertainty = jecUnc->getUncertainty(true);
652	<	jet.JEC_factor_raw = pat_jet.jecFactor("Uncorrected");
653	<
654	<	jet.btag_simpleSecondaryVertexHighEff=pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
655	<	jet.btag_simpleSecondaryVertexHighPur=pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
656	<	jet.btag_combinedSecondaryVertex=pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags");
657	<	jet.btag_combinedSecondaryVertexMVA=pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
658	<	jet.btag_jetBProbability=pat_jet.bDiscriminator("jetBProbabilityBJetTags");
659	<	jet.btag_jetProbability=pat_jet.bDiscriminator("jetProbabilityBJetTags");
636	>	jet.set_neutralEmEnergyFraction (pat_jet.neutralEmEnergyFraction());
637	>	jet.set_neutralHadronEnergyFraction (pat_jet.neutralHadronEnergyFraction());
638	>	jet.set_chargedEmEnergyFraction (pat_jet.chargedEmEnergyFraction());
639	>	jet.set_chargedHadronEnergyFraction (pat_jet.chargedHadronEnergyFraction());
640	>	jet.set_muonEnergyFraction (pat_jet.muonEnergyFraction());
641	>	jet.set_photonEnergyFraction (pat_jet.photonEnergyFraction());
642	>	jet.set_chargedMultiplicity (pat_jet.chargedMultiplicity());
643	>	jet.set_neutralMultiplicity (pat_jet.neutralMultiplicity());
644	>	jet.set_muonMultiplicity (pat_jet.muonMultiplicity());
645	>	jet.set_electronMultiplicity (pat_jet.electronMultiplicity());
646	>	jet.set_photonMultiplicity (pat_jet.photonMultiplicity());
647	>	}
648	>	if(doJECUncertainty){
649	>	jecUnc->setJetEta(pat_jet.eta());
650	>	jecUnc->setJetPt(pat_jet.pt());
651	>	jet.set_JEC_uncertainty(jecUnc->getUncertainty(true));
652	>	}
653	>	jet.set_JEC_factor_raw(pat_jet.jecFactor("Uncorrected"));
654	>
655	>	jet.set_btag_simpleSecondaryVertexHighEff(pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
656	>	jet.set_btag_simpleSecondaryVertexHighPur(pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
657	>	jet.set_btag_combinedSecondaryVertex(pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags"));
658	>	jet.set_btag_combinedSecondaryVertexMVA(pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
659	>	jet.set_btag_jetBProbability(pat_jet.bDiscriminator("jetBProbabilityBJetTags"));
660	>	jet.set_btag_jetProbability(pat_jet.bDiscriminator("jetProbabilityBJetTags"));
661
662	+
663		const reco::GenJet *genj = pat_jet.genJet();
664		if(genj){
665	<	jet.genjet_pt = genj->pt();
666	<	jet.genjet_eta = genj->eta();
667	<	jet.genjet_phi = genj->phi();
668	<	jet.genjet_energy = genj->energy();
665	>
666	>	for(unsigned int k=0; k<genjets->size(); ++k){
667	>	if(genj->pt()==genjets->at(k).pt() && genj->eta()==genjets->at(k).eta())
668	>	jet.set_genjet_index(k);
669	>	}
670	>	//         if( jet.genjet_index()<0){
671	>	//           std::cout<< "genjet not found for " << genj->pt() << "  " << genj->eta() << std::endl;
672	>	//         }
673	>
674		if(doAllGenParticles){
675		std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
676		for(unsigned int l = 0; l<jetgenps.size(); ++l){
677		for(unsigned int k=0; k< genps.size(); ++k){
678	<	if(jetgenps[l]->pt() == genps[k].pt && jetgenps[l]->pdgId() == genps[k].pdgId){
679	<	jet.genparticles_indices.push_back(genps[k].index);
678	>	if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
679	>	jet.add_genparticles_index(genps[k].index());
680	>	break;
681		}
682		}
683		}
684	<	if(jet.genparticles_indices.size()!= jetgenps.size())
685	<	std::cout << "WARNING: Found only " << jet.genparticles_indices.size() << " from " << jetgenps.size() << " gen particles of this jet"<<std::endl;
684	>	if(jet.genparticles_indices().size()!= jetgenps.size())
685	>	std::cout << "WARNING: Found only " << jet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this jet"<<std::endl;
686		}
687	+
688		}
689	+
690		jets[j].push_back(jet);
691		}
692		}
693		}
694
695	+
696	+
697		// ------------- top jets -------------
698		if(doTopJets){
699		for(size_t j=0; j< topjet_sources.size(); ++j){
#	Line 617 | Line 711 | NtupleWriter::analyze(const edm::Event&
711		if(fabs(pat_topjet.eta()) > topjet_etamax) continue;
712
713		TopJet topjet;
714	<	topjet.charge = pat_topjet.charge();
715	<	topjet.pt = pat_topjet.pt();
716	<	topjet.eta = pat_topjet.eta();
717	<	topjet.phi = pat_topjet.phi();
718	<	topjet.energy = pat_topjet.energy();
719	<	topjet.numberOfDaughters =pat_topjet.numberOfDaughters();
714	>	topjet.set_charge(pat_topjet.charge());
715	>	topjet.set_pt(pat_topjet.pt());
716	>	topjet.set_eta(pat_topjet.eta());
717	>	topjet.set_phi(pat_topjet.phi());
718	>	topjet.set_energy(pat_topjet.energy());
719	>	topjet.set_numberOfDaughters(pat_topjet.numberOfDaughters());
720		const reco::TrackRefVector&  topjettracks = pat_topjet.associatedTracks();
721	<	topjet.nTracks = topjettracks.size();
722	<	topjet.jetArea = pat_topjet.jetArea();
723	<	topjet.pileup = pat_topjet.pileup();
724	<	//       topjet.neutralEmEnergyFraction =pat_topjet.neutralEmEnergyFraction();
725	<	//       topjet.neutralHadronEnergyFraction =pat_topjet.neutralHadronEnergyFraction();
726	<	//       topjet.chargedEmEnergyFraction =pat_topjet.chargedEmEnergyFraction();
727	<	//       topjet.chargedHadronEnergyFraction =pat_topjet.chargedHadronEnergyFraction();
728	<	//       topjet.muonEnergyFraction =pat_topjet.muonEnergyFraction();
729	<	//       topjet.photonEnergyFraction =pat_topjet.photonEnergyFraction();
730	<	//       topjet.chargedMultiplicity =pat_topjet.chargedMultiplicity();
731	<	//       topjet.neutralMultiplicity =pat_topjet.neutralMultiplicity();
732	<	//       topjet.muonMultiplicity =pat_topjet.muonMultiplicity();
733	<	//       topjet.electronMultiplicity =pat_topjet.electronMultiplicity();
734	<	//       topjet.photonMultiplicity =pat_topjet.photonMultiplicity();
735	<
736	<	jecUnc->setJetEta(pat_topjet.eta());
737	<	jecUnc->setJetPt(pat_topjet.pt());
738	<	topjet.JEC_uncertainty = jecUnc->getUncertainty(true);
739	<	topjet.JEC_factor_raw = pat_topjet.jecFactor("Uncorrected");
740	<
647	<	topjet.btag_simpleSecondaryVertexHighEff=pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
648	<	topjet.btag_simpleSecondaryVertexHighPur=pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
649	<	topjet.btag_combinedSecondaryVertex=pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags");
650	<	topjet.btag_combinedSecondaryVertexMVA=pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
651	<	topjet.btag_jetBProbability=pat_topjet.bDiscriminator("jetBProbabilityBJetTags");
652	<	topjet.btag_jetProbability=pat_topjet.bDiscriminator("jetProbabilityBJetTags");
721	>	topjet.set_nTracks( topjettracks.size());
722	>	topjet.set_jetArea( pat_topjet.jetArea());
723	>	topjet.set_pileup( pat_topjet.pileup());
724	>	//       topjet.set_neutralEmEnergyFraction(pat_topjet.neutralEmEnergyFraction());
725	>	//       topjet.set_neutralHadronEnergyFraction(pat_topjet.neutralHadronEnergyFraction());
726	>	//       topjet.set_chargedEmEnergyFraction(pat_topjet.chargedEmEnergyFraction());
727	>	//       topjet.set_chargedHadronEnergyFraction(pat_topjet.chargedHadronEnergyFraction());
728	>	//       topjet.set_muonEnergyFraction(pat_topjet.muonEnergyFraction());
729	>	//       topjet.set_photonEnergyFraction(pat_topjet.photonEnergyFraction());
730	>	//       topjet.set_chargedMultiplicity(pat_topjet.chargedMultiplicity());
731	>	//       topjet.set_neutralMultiplicity(pat_topjet.neutralMultiplicity());
732	>	//       topjet.set_muonMultiplicity(pat_topjet.muonMultiplicity());
733	>	//       topjet.set_electronMultiplicity(pat_topjet.electronMultiplicity());
734	>	//       topjet.set_photonMultiplicity(pat_topjet.photonMultiplicity());
735	>	if(doJECUncertainty){
736	>	jecUnc->setJetEta(pat_topjet.eta());
737	>	jecUnc->setJetPt(pat_topjet.pt());
738	>	topjet.set_JEC_uncertainty( jecUnc->getUncertainty(true));
739	>	}
740	>	topjet.set_JEC_factor_raw( pat_topjet.jecFactor("Uncorrected"));
741
742	+	topjet.set_btag_simpleSecondaryVertexHighEff(pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
743	+	topjet.set_btag_simpleSecondaryVertexHighPur(pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
744	+	topjet.set_btag_combinedSecondaryVertex(pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags"));
745	+	topjet.set_btag_combinedSecondaryVertexMVA(pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
746	+	topjet.set_btag_jetBProbability(pat_topjet.bDiscriminator("jetBProbabilityBJetTags"));
747	+	topjet.set_btag_jetProbability(pat_topjet.bDiscriminator("jetProbabilityBJetTags"));
748	+
749	+	/*
750		const reco::GenJet *genj = pat_topjet.genJet();
751		if(genj){
752	<	topjet.genjet_pt = genj->pt();
753	<	topjet.genjet_eta = genj->eta();
754	<	topjet.genjet_phi = genj->phi();
755	<	topjet.genjet_energy = genj->energy();
752	>	topjet.set_genjet_pt ( genj->pt());
753	>	topjet.set_genjet_eta ( genj->eta());
754	>	topjet.set_genjet_phi ( genj->phi());
755	>	topjet.set_genjet_energy ( genj->energy());
756		if(doAllGenParticles){
757		std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
758		for(unsigned int l = 0; l<jetgenps.size(); ++l){
759		for(unsigned int k=0; k< genps.size(); ++k){
760	<	if(jetgenps[l]->pt() == genps[k].pt && jetgenps[l]->pdgId() == genps[k].pdgId){
761	<	topjet.genparticles_indices.push_back(genps[k].index);
760	>	if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
761	>	topjet.add_genparticles_index(genps[k].index());
762	>	break;
763		}
764		}
765		}
766	<	if(topjet.genparticles_indices.size()!= jetgenps.size())
767	<	std::cout << "WARNING: Found only " << topjet.genparticles_indices.size() << " from " << jetgenps.size() << " gen particles of this topjet"<<std::endl;
766	>	if(topjet.genparticles_indices().size()!= jetgenps.size())
767	>	std::cout << "WARNING: Found only " << topjet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this topjet"<<std::endl;
768		}
769		}
770	+	*/
771
772		for (unsigned int k = 0; k < pat_topjet.numberOfDaughters(); k++) {
773		Particle subjet_v4;
774	<	subjet_v4.pt = pat_topjet.daughter(k)->p4().pt();
775	<	subjet_v4.eta = pat_topjet.daughter(k)->p4().eta();
776	<	subjet_v4.phi = pat_topjet.daughter(k)->p4().phi();
777	<	subjet_v4.energy = pat_topjet.daughter(k)->p4().E();
778	<	topjet.subjets.push_back(subjet_v4);
774	>
775	>	reco::Candidate const * subjetd =  pat_topjet.daughter(k);
776	>	pat::Jet const * patsubjetd = dynamic_cast<pat::Jet const *>(subjetd);
777	>	if(patsubjetd)
778	>	{
779	>	subjet_v4.set_pt(patsubjetd->correctedP4(0).pt());
780	>	subjet_v4.set_eta(patsubjetd->correctedP4(0).eta());
781	>	subjet_v4.set_phi(patsubjetd->correctedP4(0).phi());
782	>	subjet_v4.set_energy(patsubjetd->correctedP4(0).E());
783	>	topjet.add_subjet(subjet_v4);
784	>	//btag info
785	>	topjet.add_subFlavour(patsubjetd->partonFlavour());
786	>	topjet.add_subCSV(patsubjetd->bDiscriminator("combinedSecondaryVertexBJetTags"));
787	>	if (doTagInfos)
788	>	{
789	>	//ip tag info
790	>	reco::TaggingVariableList tvlIP=patsubjetd->tagInfoTrackIP("impactParameter")->taggingVariables();
791	>	topjet.add_subTrackMomentum(tvlIP.getList(reco::btau::trackMomentum,false));
792	>	topjet.add_subTrackEta(tvlIP.getList(reco::btau::trackEta,false));
793	>	topjet.add_subTrackEtaRel(tvlIP.getList(reco::btau::trackEtaRel,false));
794	>	topjet.add_subTrackDeltaR(tvlIP.getList(reco::btau::trackDeltaR,false));
795	>	topjet.add_subTrackSip3dVal(tvlIP.getList(reco::btau::trackSip3dVal,false));
796	>	topjet.add_subTrackSip3dSig(tvlIP.getList(reco::btau::trackSip3dSig,false));
797	>	topjet.add_subTrackSip2dVal(tvlIP.getList(reco::btau::trackSip2dVal,false));
798	>	topjet.add_subTrackSip2dSig(tvlIP.getList(reco::btau::trackSip2dSig,false));
799	>	topjet.add_subTrackDecayLenVal(tvlIP.getList(reco::btau::trackDecayLenVal,false));
800	>	topjet.add_subTrackChi2(tvlIP.getList(reco::btau::trackChi2,false));
801	>	topjet.add_subTrackNTotalHits(tvlIP.getList(reco::btau::trackNTotalHits,false));
802	>	topjet.add_subTrackNPixelHits(tvlIP.getList(reco::btau::trackNPixelHits,false));
803	>	topjet.add_subTrackPtRel(tvlIP.getList(reco::btau::trackPtRel,false));
804	>	topjet.add_subTrackPPar(tvlIP.getList(reco::btau::trackPPar,false));
805	>	topjet.add_subTrackPtRatio(tvlIP.getList(reco::btau::trackPtRatio,false));
806	>	topjet.add_subTrackPParRatio(tvlIP.getList(reco::btau::trackPParRatio,false));
807	>	topjet.add_subTrackJetDistVal(tvlIP.getList(reco::btau::trackJetDistVal,false));
808	>	topjet.add_subTrackJetDistSig(tvlIP.getList(reco::btau::trackJetDistSig,false));
809	>	topjet.add_subTrackGhostTrackDistVal(tvlIP.getList(reco::btau::trackGhostTrackDistVal,false));
810	>	topjet.add_subTrackGhostTrackDistSig(tvlIP.getList(reco::btau::trackGhostTrackDistSig,false));
811	>	topjet.add_subTrackGhostTrackWeight(tvlIP.getList(reco::btau::trackGhostTrackWeight,false));
812	>	//sv tag info
813	>	reco::TaggingVariableList tvlSV=patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->taggingVariables();
814	>	topjet.add_subFlightDistance2dVal(tvlSV.getList(reco::btau::flightDistance2dVal,false));
815	>	topjet.add_subFlightDistance2dSig(tvlSV.getList(reco::btau::flightDistance2dSig,false));
816	>	topjet.add_subFlightDistance3dVal(tvlSV.getList(reco::btau::flightDistance3dVal,false));
817	>	topjet.add_subFlightDistance3dSig(tvlSV.getList(reco::btau::flightDistance3dSig,false));
818	>	topjet.add_subVertexJetDeltaR(tvlSV.getList(reco::btau::vertexJetDeltaR,false));
819	>	topjet.add_subJetNSecondaryVertices(tvlSV.get(reco::btau::jetNSecondaryVertices,-9999));
820	>	topjet.add_subVertexNTracks(tvlSV.get(reco::btau::vertexNTracks,-9999));
821	>	std::vector<TLorentzVector> vp4; vp4.clear();
822	>	std::vector<float> vchi2; vchi2.clear();
823	>	std::vector<float> vndof; vndof.clear();
824	>	std::vector<float> vchi2ndof; vchi2ndof.clear();
825	>	std::vector<float> vtrsize; vtrsize.clear();
826	>	for(unsigned int i=0; i<patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->nVertices(); i++)
827	>	{
828	>	ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > p4 = patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).p4();
829	>	vp4.push_back(TLorentzVector(p4.px(),p4.py(),p4.pz(),p4.e()));
830	>	vchi2.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).chi2());
831	>	vndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).ndof());
832	>	vchi2ndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).normalizedChi2());
833	>	vtrsize.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).tracksSize());
834	>	}
835	>	topjet.add_subSecondaryVertex(vp4);
836	>	topjet.add_subVertexChi2(vchi2);
837	>	topjet.add_subVertexNdof(vndof);
838	>	topjet.add_subVertexNormalizedChi2(vchi2ndof);
839	>	topjet.add_subVertexTracksSize(vtrsize);
840	>	//try computer
841	>	edm::ESHandle<JetTagComputer> computerHandle;
842	>	iSetup.get<JetTagComputerRecord>().get( SVComputer_.label(), computerHandle );
843	>	const GenericMVAJetTagComputer *computer = dynamic_cast<const GenericMVAJetTagComputer*>( computerHandle.product() );
844	>	if(computer)
845	>	{
846	>	computer->passEventSetup(iSetup);
847	>	std::vector<const reco::BaseTagInfo*>  baseTagInfos;
848	>	baseTagInfos.push_back(patsubjetd->tagInfoTrackIP("impactParameter") );
849	>	baseTagInfos.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex") );
850	>	JetTagComputer::TagInfoHelper helper(baseTagInfos);
851	>	reco::TaggingVariableList vars = computer->taggingVariables(helper);
852	>	topjet.add_subVertexMassJTC(vars.get(reco::btau::vertexMass,-9999));
853	>	topjet.add_subVertexCategoryJTC(vars.get(reco::btau::vertexCategory,-9999));
854	>	topjet.add_subVertexEnergyRatioJTC(vars.get(reco::btau::vertexEnergyRatio,-9999));
855	>	topjet.add_subTrackSip3dSigAboveCharmJTC(vars.get(reco::btau::trackSip3dSigAboveCharm,-9999));
856	>	}
857	>	}
858	>	}
859	>	else
860	>	{
861	>	//filling only standard information in case the subjet has not been pat-tified during the pattuples production
862	>	subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt());
863	>	subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta());
864	>	subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi());
865	>	subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E());
866	>	topjet.add_subjet(subjet_v4);
867	>	}
868	>
869	>
870		}
871		topjets[j].push_back(topjet);
872		}
873		}
874		}
875	<
876	<
875	>
876	>
877		// ------------- generator top jets -------------
878		if(doGenTopJets){
879		for(size_t j=0; j< gentopjet_sources.size(); ++j){
880
881		gentopjets[j].clear();
882	<
882	>
883		edm::Handle<reco::BasicJetCollection> reco_gentopjets;
884		//edm::Handle<std::vector<reco::Jet> > reco_gentopjets;
885		iEvent.getByLabel(gentopjet_sources[j], reco_gentopjets);
#	Line 702 | Line 891 | NtupleWriter::analyze(const edm::Event&
891		if(fabs(reco_gentopjet.eta()) > gentopjet_etamax) continue;
892
893		TopJet gentopjet;
894	<	gentopjet.charge = reco_gentopjet.charge();
895	<	gentopjet.pt = reco_gentopjet.pt();
896	<	gentopjet.eta = reco_gentopjet.eta();
897	<	gentopjet.phi = reco_gentopjet.phi();
898	<	gentopjet.energy = reco_gentopjet.energy();
899	<	gentopjet.numberOfDaughters =reco_gentopjet.numberOfDaughters();
894	>	gentopjet.set_charge(reco_gentopjet.charge());
895	>	gentopjet.set_pt(reco_gentopjet.pt());
896	>	gentopjet.set_eta(reco_gentopjet.eta());
897	>	gentopjet.set_phi(reco_gentopjet.phi());
898	>	gentopjet.set_energy(reco_gentopjet.energy());
899	>	gentopjet.set_numberOfDaughters(reco_gentopjet.numberOfDaughters());
900
901		for (unsigned int k = 0; k < reco_gentopjet.numberOfDaughters(); k++) {
902		Particle subjet_v4;
903	<	subjet_v4.pt = reco_gentopjet.daughter(k)->p4().pt();
904	<	subjet_v4.eta = reco_gentopjet.daughter(k)->p4().eta();
905	<	subjet_v4.phi = reco_gentopjet.daughter(k)->p4().phi();
906	<	subjet_v4.energy = reco_gentopjet.daughter(k)->p4().E();
907	<	gentopjet.subjets.push_back(subjet_v4);
903	>	subjet_v4.set_pt(reco_gentopjet.daughter(k)->p4().pt());
904	>	subjet_v4.set_eta(reco_gentopjet.daughter(k)->p4().eta());
905	>	subjet_v4.set_phi(reco_gentopjet.daughter(k)->p4().phi());
906	>	subjet_v4.set_energy(reco_gentopjet.daughter(k)->p4().E());
907	>	gentopjet.add_subjet(subjet_v4);
908		}
909		gentopjets[j].push_back(gentopjet);
910		}
#	Line 734 | Line 923 | NtupleWriter::analyze(const edm::Event&
923		for (unsigned int i = 0; i < pat_photons.size(); ++i) {
924		pat::Photon pat_photon = pat_photons[i];
925		Photon ph;
926	<	ph.charge = 0;
927	<	ph.pt =  pat_photon.pt();
928	<	ph.eta =  pat_photon.eta();
929	<	ph.phi =  pat_photon.phi();
930	<	ph.energy =  pat_photon.energy();
931	<	ph.vertex_x = pat_photon.vertex().x();
932	<	ph.vertex_y = pat_photon.vertex().y();
933	<	ph.vertex_z = pat_photon.vertex().z();
934	<	ph.supercluster_eta = pat_photon.superCluster()->eta();
935	<	ph.supercluster_phi = pat_photon.superCluster()->phi();
936	<	//       ph.neutralHadronIso = pat_photon.neutralHadronIso();
937	<	//       ph.chargedHadronIso = pat_photon.chargedHadronIso();
938	<	ph.trackIso = pat_photon.trackIso();
926	>	ph.set_charge(0);
927	>	ph.set_pt( pat_photon.pt());
928	>	ph.set_eta( pat_photon.eta());
929	>	ph.set_phi( pat_photon.phi());
930	>	ph.set_energy( pat_photon.energy());
931	>	ph.set_vertex_x(pat_photon.vertex().x());
932	>	ph.set_vertex_y(pat_photon.vertex().y());
933	>	ph.set_vertex_z(pat_photon.vertex().z());
934	>	ph.set_supercluster_eta(pat_photon.superCluster()->eta());
935	>	ph.set_supercluster_phi(pat_photon.superCluster()->phi());
936	>	//       ph.set_neutralHadronIso(pat_photon.neutralHadronIso());
937	>	//       ph.set_chargedHadronIso(pat_photon.chargedHadronIso());
938	>	ph.set_trackIso(pat_photon.trackIso());
939		phs[j].push_back(ph);
940		}
941		}
#	Line 766 | Line 955 | NtupleWriter::analyze(const edm::Event&
955		else{
956		pat::MET pat_met = pat_mets[0];
957
958	<	met[j].pt=pat_met.pt();
959	<	met[j].phi=pat_met.phi();
960	<	met[j].mEtSig= pat_met.mEtSig();
958	>	met[j].set_pt(pat_met.pt());
959	>	met[j].set_phi(pat_met.phi());
960	>	met[j].set_mEtSig(pat_met.mEtSig());
961		}
962
963		}
#	Line 854 | Line 1043 | void
1043		NtupleWriter::beginRun(edm::Run const& iRun, edm::EventSetup const&  iSetup)
1044		{
1045		if(doTrigger){
1046	<	bool setup_changed = false;
1047	<	hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
1046	>	//bool setup_changed = false;
1047	>	//hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
1048		newrun=true;
1049		}
1050
1051		if(doJets || doTopJets){
1052	<	edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
1053	<	iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
1054	<	JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
1055	<	jecUnc = new JetCorrectionUncertainty(JetCorPar);
1052	>	if(doJECUncertainty){
1053	>	edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
1054	>	iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
1055	>	JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
1056	>	jecUnc = new JetCorrectionUncertainty(JetCorPar);
1057	>	}
1058		}
1059		}
1060

Diff Legend

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