[ViewVC] Diff of: cvsroot/UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc

Comparing UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc (file contents):
Revision 1.6 by peiffer, Wed Apr 11 15:15:44 2012 UTC vs.
Revision 1.26 by peiffer, Mon Jun 10 13:33:05 2013 UTC

+//
+//
-<
+#include "UHHAnalysis//NtupleWriter/interface/NtupleWriter.h"
-<
-<
->
+#include "UHHAnalysis/NtupleWriter/interface/NtupleWriter.h"
->
+#include "RecoBTau/JetTagComputer/interface/GenericMVAJetTagComputer.h"
->
+#include "RecoBTau/JetTagComputer/interface/GenericMVAJetTagComputerWrapper.h"
->
+#include "RecoBTau/JetTagComputer/interface/JetTagComputer.h"
->
+#include "RecoBTau/JetTagComputer/interface/JetTagComputerRecord.h"
->
+#include "RecoBTag/SecondaryVertex/interface/CombinedSVComputer.h"
+//
+// constants, enums and typedefs
+  doMuons = iConfig.getParameter<bool>("doMuons");
+  doTaus = iConfig.getParameter<bool>("doTaus");
+  doJets = iConfig.getParameter<bool>("doJets");
-+
+  doGenJets = iConfig.getParameter<bool>("doGenJets");
-+
+  doJECUncertainty = iConfig.getParameter<bool>("doJECUncertainty");
-+
+  doGenTopJets = iConfig.getParameter<bool>("doGenTopJets");
+  doPhotons = iConfig.getParameter<bool>("doPhotons");
+  doMET = iConfig.getParameter<bool>("doMET");
+  doGenInfo = iConfig.getParameter<bool>("doGenInfo");
-+
+  doAllGenParticles = iConfig.getParameter<bool>("doAllGenParticles");
-+
+  doLumiInfo = iConfig.getParameter<bool>("doLumiInfo");
+  doPV = iConfig.getParameter<bool>("doPV");
+  doTopJets = iConfig.getParameter<bool>("doTopJets");
+  doTrigger = iConfig.getParameter<bool>("doTrigger");
-<
->
+  SVComputer_  = iConfig.getUntrackedParameter<edm::InputTag>("svComputer",edm::InputTag("combinedSecondaryVertex"));
->
+  doTagInfos = iConfig.getUntrackedParameter<bool>("doTagInfos",false);
+  // initialization of tree variables
+  tr->Branch("run",&run);
+  tr->Branch("event",&event);
+  tr->Branch("luminosityBlock",&luminosityBlock);
+  tr->Branch("isRealData",&isRealData);
-<
+  tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
-<
+  tr->Branch("beamspot_x0",&beamspot_x0);
-<
+  tr->Branch("beamspot_y0",&beamspot_y0);
-<
+  tr->Branch("beamspot_z0",&beamspot_z0);
->
+  tr->Branch("rho",&rho);
->
+  rho_source = iConfig.getParameter<edm::InputTag>("rho_source");
-+
+  //tr->Branch("HBHENoiseFilterResult",&HBHENoiseFilterResult);
-+
+  if(doLumiInfo){
-+
+    tr->Branch("intgRecLumi",&intgRecLumi);
-+
+    tr->Branch("intgDelLumi",&intgDelLumi);
-+
+  }
-+
+  if(doPV){
-+
+    tr->Branch("beamspot_x0",&beamspot_x0);
-+
+    tr->Branch("beamspot_y0",&beamspot_y0);
-+
+    tr->Branch("beamspot_z0",&beamspot_z0);
-+
+  }
+  if(doElectrons){
+    electron_sources = iConfig.getParameter<std::vector<std::string> >("electron_sources");
+    for(size_t j=0; j< electron_sources.size(); ++j){
+      tr->Branch( jet_sources[j].c_str(), "std::vector<Jet>", &jets[j]);
-+
+  if(doGenJets){
-+
+    genjet_sources = iConfig.getParameter<std::vector<std::string> >("genjet_sources");
-+
+    genjet_ptmin = iConfig.getParameter<double> ("genjet_ptmin");
-+
+    genjet_etamax = iConfig.getParameter<double> ("genjet_etamax");
-+
+    for(size_t j=0; j< genjet_sources.size(); ++j){
-+
+      tr->Branch( genjet_sources[j].c_str(), "std::vector<Particle>", &genjets[j]);
-+
+    }
-+
+  }
+  if(doTopJets){
+    topjet_sources = iConfig.getParameter<std::vector<std::string> >("topjet_sources");
+    topjet_ptmin = iConfig.getParameter<double> ("topjet_ptmin");
+      tr->Branch( topjet_sources[j].c_str(), "std::vector<TopJet>", &topjets[j]);
-+
+  if(doGenTopJets){
-+
+    gentopjet_sources = iConfig.getParameter<std::vector<std::string> >("gentopjet_sources");
-+
+    gentopjet_ptmin = iConfig.getParameter<double> ("gentopjet_ptmin");
-+
+    gentopjet_etamax = iConfig.getParameter<double> ("gentopjet_etamax");
-+
+    for(size_t j=0; j< gentopjet_sources.size(); ++j){
-+
+      tr->Branch( gentopjet_sources[j].c_str(), "std::vector<TopJet>", &gentopjets[j]);
-+
+    }
-+
+  }
+  if(doPhotons){
+    photon_sources = iConfig.getParameter<std::vector<std::string> >("photon_sources");
+    for(size_t j=0; j< photon_sources.size(); ++j){
+  if(doGenInfo){
-+
+    genparticle_source= iConfig.getParameter<edm::InputTag>("genparticle_source");
+    tr->Branch("genInfo","GenInfo",&genInfo);
+    tr->Branch("GenParticles","std::vector<GenParticle>", &genps);
+    //tr->Branch("triggerResults","std::map<std::string, bool>",&triggerResults);
+    tr->Branch("triggerNames", "std::vector<std::string>", &triggerNames);
+    tr->Branch("triggerResults", "std::vector<bool>", &triggerResults);
-<
+    tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
-<
+    tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
->
+    //tr->Branch("L1_prescale", "std::vector<int>", &L1_prescale);
->
+    //tr->Branch("HLT_prescale", "std::vector<int>", &HLT_prescale);
+  newrun = true;
+   luminosityBlock = iEvent.luminosityBlock();
+   isRealData      = iEvent.isRealData();
-<
+   if(isRealData){
-<
+     edm::Handle<bool> bool_handle;
-<
+     iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
-<
+     HBHENoiseFilterResult = *(bool_handle.product());
-<
+   }
-<
+   else HBHENoiseFilterResult = false;
->
+   edm::Handle<double> m_rho;
->
+   iEvent.getByLabel(rho_source,m_rho);
->
+   rho=*m_rho;
->
->
+//    if(isRealData){
->
+//      edm::Handle<bool> bool_handle;
->
+//      iEvent.getByLabel(edm::InputTag("HBHENoiseFilterResultProducer","HBHENoiseFilterResult"),bool_handle);
->
+//      HBHENoiseFilterResult = *(bool_handle.product());
->
+//    }
->
+//    else HBHENoiseFilterResult = false;
+   // ------------- primary vertices and beamspot  -------------
+         reco::Vertex reco_pv = reco_pvs[i];
+         PrimaryVertex pv;
-<
+         pv.x =  reco_pv.x();
-<
+         pv.y =  reco_pv.y();
-<
+         pv.z =  reco_pv.z();
-<
+         pv.nTracks =  reco_pv.nTracks();
-<
+         //pv.isValid =  reco_pv.isValid();
-<
+         pv.chi2 =  reco_pv.chi2();
-<
+         pv.ndof =  reco_pv.ndof();
->
+         pv.set_x( reco_pv.x());
->
+         pv.set_y( reco_pv.y());
->
+         pv.set_z( reco_pv.z());
->
+         pv.set_nTracks( reco_pv.nTracks());
->
+         //pv.set_isValid( reco_pv.isValid());
->
+         pv.set_chi2( reco_pv.chi2());
->
+         pv.set_ndof( reco_pv.ndof());
+         pvs[j].push_back(pv);
-+
-+
+     edm::Handle<reco::BeamSpot> beamSpot;
-+
+     iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
-+
+     const reco::BeamSpot & bsp = *beamSpot;
-+
-+
+     beamspot_x0 = bsp.x0();
-+
+     beamspot_y0 = bsp.y0();
-+
+     beamspot_z0 = bsp.z0();
-+
-+
+ // ------------- generator info -------------
-<
+   edm::Handle<reco::BeamSpot> beamSpot;
-<
+   iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
-<
+   const reco::BeamSpot & bsp = *beamSpot;
-<
-<
+   beamspot_x0 = bsp.x0();
-<
+   beamspot_y0 = bsp.y0();
-<
+   beamspot_z0 = bsp.z0();
->
+   if(doGenInfo){
->
+     genInfo.clear_weights();
->
+     genInfo.clear_binningValues();
->
+     genps.clear();
->
->
+     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
->
+     iEvent.getByLabel("generator", genEventInfoProduct);
->
+     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
->
->
+     for(unsigned int k=0; k<genEventInfo.binningValues().size();++k){
->
+       genInfo.add_binningValue(genEventInfo.binningValues().at(k));
->
+     }
->
+     for(unsigned int k=0; k<genEventInfo.weights().size();++k){
->
+       genInfo.add_weight(genEventInfo.weights().at(k));
->
+     }
->
+     genInfo.set_alphaQCD(genEventInfo.alphaQCD());
->
+     genInfo.set_alphaQED(genEventInfo.alphaQED());
->
+     genInfo.set_qScale(genEventInfo.qScale());
->
->
+     const gen::PdfInfo* pdf = genEventInfo.pdf();
->
+     if(pdf){
->
+       genInfo.set_pdf_id1(pdf->id.first);
->
+       genInfo.set_pdf_id2(pdf->id.second);
->
+       genInfo.set_pdf_x1(pdf->x.first);
->
+       genInfo.set_pdf_x2(pdf->x.second);
->
+       genInfo.set_pdf_scalePDF(pdf->scalePDF);
->
+       genInfo.set_pdf_xPDF1(pdf->xPDF.first);
->
+       genInfo.set_pdf_xPDF2(pdf->xPDF.second);
->
+     }
->
+     else{
->
+       genInfo.set_pdf_id1(-999);
->
+       genInfo.set_pdf_id2(-999);
->
+       genInfo.set_pdf_x1(-999);
->
+       genInfo.set_pdf_x2(-999);
->
+       genInfo.set_pdf_scalePDF(-999);
->
+       genInfo.set_pdf_xPDF1(-999);
->
+       genInfo.set_pdf_xPDF2(-999);
->
+     }
->
->
+     edm::Handle<std::vector<PileupSummaryInfo> > pus;
->
+     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
->
+     genInfo.set_pileup_NumInteractions_intime(0);
->
+     genInfo.set_pileup_NumInteractions_ootbefore(0);
->
+     genInfo.set_pileup_NumInteractions_ootafter(0);
->
+     if(pus.isValid()){
->
+       genInfo.set_pileup_TrueNumInteractions ( (float) pus->at(0).getTrueNumInteractions());
->
+       for(size_t i=0; i<pus->size(); ++i){
->
+         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
->
+           genInfo.set_pileup_NumInteractions_intime( genInfo.pileup_NumInteractions_intime() + pus->at(i).getPU_NumInteractions());
->
+         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
->
+           genInfo.set_pileup_NumInteractions_ootbefore( genInfo.pileup_NumInteractions_ootbefore() + pus->at(i).getPU_NumInteractions());
->
+         }
->
+         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
->
+           genInfo.set_pileup_NumInteractions_ootafter( genInfo.pileup_NumInteractions_ootafter() + pus->at(i).getPU_NumInteractions());
->
+         }
->
+       }
->
+     }
->
->
+     edm::Handle<reco::GenParticleCollection> genPartColl;
->
+     iEvent.getByLabel(genparticle_source, genPartColl);
->
+     int index=-1;
->
+     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
->
+       index++;
->
->
+       //write out only top quarks,final state leptons and status 3 particles (works fine only for MadGraph)
->
+       bool islepton = iter->status()==1 && abs(iter->pdgId())>=11 && abs(iter->pdgId())<=16 ;
->
+       if(abs(iter->pdgId())==6 || iter->status()==3 || islepton ||  doAllGenParticles){
->
+         GenParticle genp;
->
+         genp.set_charge(iter->charge());
->
+         genp.set_pt(iter->p4().pt());
->
+         genp.set_eta(iter->p4().eta());
->
+         genp.set_phi(iter->p4().phi());
->
+         genp.set_energy(iter->p4().E());
->
+         genp.set_index(index);
->
+         genp.set_status( iter->status());
->
+         genp.set_pdgId( iter->pdgId());
->
->
+         genp.set_mother1(-1);
->
+         genp.set_mother2(-1);
->
+         genp.set_daughter1(-1);
->
+         genp.set_daughter2(-1);
->
->
+         int nm=iter->numberOfMothers();
->
+         int nd=iter->numberOfDaughters();
->
->
->
+         if (nm>0) genp.set_mother1( iter->motherRef(0).key());
->
+         if (nm>1) genp.set_mother2( iter->motherRef(1).key());
->
+         if (nd>0) genp.set_daughter1( iter->daughterRef(0).key());
->
+         if (nd>1) genp.set_daughter2( iter->daughterRef(1).key());
->
->
+         genps.push_back(genp);
->
+       }
->
+     }
->
+   }
+   // ------------- electrons -------------
+   if(doElectrons){
-+
-+
+//      edm::Handle<reco::ConversionCollection> hConversions;
-+
+//      iEvent.getByLabel("allConversions", hConversions);
-+
-+
+//      edm::Handle<reco::BeamSpot> beamSpot;
-+
+//      iEvent.getByLabel(edm::InputTag("offlineBeamSpot"), beamSpot);
-+
+//      const reco::BeamSpot & bsp = *beamSpot;
-+
+     for(size_t j=0; j< electron_sources.size(); ++j){
+       eles[j].clear();
+       edm::Handle< std::vector<pat::Electron> > ele_handle;
+         pat::Electron pat_ele = pat_electrons[i];
+         Electron ele;
-<
+         ele.charge =  pat_ele.charge();
-<
+         ele.pt =  pat_ele.pt();
-<
+         ele.eta =  pat_ele.eta();
-<
+         ele.phi =  pat_ele.phi();
-<
+         ele.energy =  pat_ele.energy();
-<
+         ele.vertex_x = pat_ele.vertex().x();
-<
+         ele.vertex_y = pat_ele.vertex().y();
-<
+         ele.vertex_z = pat_ele.vertex().z();
-<
+         ele.supercluster_eta = pat_ele.superCluster()->eta();
-<
+         ele.supercluster_phi = pat_ele.superCluster()->phi();
-<
+         ele.dB = pat_ele.dB();
-<
+         //ele.particleIso = pat_ele.particleIso();
-<
+         ele.neutralHadronIso = pat_ele.neutralHadronIso();
-<
+         ele.chargedHadronIso = pat_ele.chargedHadronIso();
-<
+         ele.trackIso = pat_ele.trackIso();
-<
+         ele.puChargedHadronIso = pat_ele.puChargedHadronIso();
-<
+         ele.gsfTrack_trackerExpectedHitsInner_numberOfLostHits = pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits();
-<
+         ele.gsfTrack_px= pat_ele.gsfTrack()->px();
-<
+         ele.gsfTrack_py= pat_ele.gsfTrack()->py();
-<
+         ele.gsfTrack_pz= pat_ele.gsfTrack()->pz();
-<
+         ele.gsfTrack_vx= pat_ele.gsfTrack()->vx();
-<
+         ele.gsfTrack_vy= pat_ele.gsfTrack()->vy();
-<
+         ele.gsfTrack_vz= pat_ele.gsfTrack()->vz();
->
+         ele.set_charge( pat_ele.charge());
->
+         ele.set_pt( pat_ele.pt());
->
+         ele.set_eta( pat_ele.eta());
->
+         ele.set_phi( pat_ele.phi());
->
+         ele.set_energy( pat_ele.energy());
->
+         ele.set_vertex_x(pat_ele.vertex().x());
->
+         ele.set_vertex_y(pat_ele.vertex().y());
->
+         ele.set_vertex_z(pat_ele.vertex().z());
->
+         ele.set_supercluster_eta(pat_ele.superCluster()->eta());
->
+         ele.set_supercluster_phi(pat_ele.superCluster()->phi());
->
+         ele.set_dB(pat_ele.dB());
->
+         //ele.set_particleIso(pat_ele.particleIso());
->
+         ele.set_neutralHadronIso(pat_ele.neutralHadronIso());
->
+         ele.set_chargedHadronIso(pat_ele.chargedHadronIso());
->
+         ele.set_trackIso(pat_ele.trackIso());
->
+         ele.set_photonIso(pat_ele.photonIso());
->
+         ele.set_puChargedHadronIso(pat_ele.puChargedHadronIso());
->
+         ele.set_gsfTrack_trackerExpectedHitsInner_numberOfLostHits(pat_ele.gsfTrack()->trackerExpectedHitsInner().numberOfLostHits());
->
+         ele.set_gsfTrack_px( pat_ele.gsfTrack()->px());
->
+         ele.set_gsfTrack_py( pat_ele.gsfTrack()->py());
->
+         ele.set_gsfTrack_pz( pat_ele.gsfTrack()->pz());
->
+         ele.set_gsfTrack_vx( pat_ele.gsfTrack()->vx());
->
+         ele.set_gsfTrack_vy( pat_ele.gsfTrack()->vy());
->
+         ele.set_gsfTrack_vz( pat_ele.gsfTrack()->vz());
->
+         //ele.set_passconversionveto(!ConversionTools::hasMatchedConversion(pat_ele,hConversions,bsp.position()));
->
+         ele.set_passconversionveto(pat_ele.passConversionVeto());
->
+         ele.set_dEtaIn(pat_ele.deltaEtaSuperClusterTrackAtVtx());
->
+         ele.set_dPhiIn(pat_ele.deltaPhiSuperClusterTrackAtVtx());
->
+         ele.set_sigmaIEtaIEta(pat_ele.sigmaIetaIeta());
->
+         ele.set_HoverE(pat_ele.hadronicOverEm());
->
+         ele.set_fbrem(pat_ele.fbrem());
->
+         ele.set_EoverPIn(pat_ele.eSuperClusterOverP());
->
+         ele.set_EcalEnergy(pat_ele.ecalEnergy());
->
+         ele.set_mvaTrigV0(pat_ele.electronID("mvaTrigV0"));
->
+         ele.set_mvaNonTrigV0(pat_ele.electronID("mvaNonTrigV0"));
->
+         float AEff03 = 0.00;
->
+         if(isRealData){
->
+           AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAData2011);
->
+         }else{
->
+           AEff03 = ElectronEffectiveArea::GetElectronEffectiveArea(ElectronEffectiveArea::kEleGammaAndNeutralHadronIso03, pat_ele.superCluster()->eta(), ElectronEffectiveArea::kEleEAFall11MC);
->
+         }
->
+         ele.set_AEff(AEff03);
->
+         eles[j].push_back(ele);
+         pat::Muon pat_mu = pat_muons[i];
+         Muon mu;
-<
+         mu.charge =  pat_mu.charge();
-<
+         mu.pt =  pat_mu.pt();
-<
+         mu.eta =  pat_mu.eta();
-<
+         mu.phi =  pat_mu.phi();
-<
+         mu.energy =  pat_mu.energy();
-<
+         mu.vertex_x = pat_mu.vertex().x();
-<
+         mu.vertex_y = pat_mu.vertex().y();
-<
+         mu.vertex_z = pat_mu.vertex().z();
-<
+         mu.dB = pat_mu.dB();
-<
+         //mu.particleIso = pat_mu.particleIso();
-<
+         mu.neutralHadronIso = pat_mu.neutralHadronIso();
-<
+         mu.chargedHadronIso = pat_mu.chargedHadronIso();
-<
+         mu.trackIso = pat_mu.trackIso();
-<
+         mu.puChargedHadronIso = pat_mu.puChargedHadronIso();
-<
+         mu.isGlobalMuon = pat_mu.isGlobalMuon();
-<
+         mu.isStandAloneMuon = pat_mu.isStandAloneMuon();
-<
+         mu.isTrackerMuon = pat_mu.isTrackerMuon();
-<
+         mu.numberOfMatchedStations = pat_mu.numberOfMatchedStations();
->
+         mu.set_charge( pat_mu.charge());
->
+         mu.set_pt( pat_mu.pt());
->
+         mu.set_eta( pat_mu.eta());
->
+         mu.set_phi( pat_mu.phi());
->
+         mu.set_energy( pat_mu.energy());
->
+         mu.set_vertex_x ( pat_mu.vertex().x());
->
+         mu.set_vertex_y ( pat_mu.vertex().y());
->
+         mu.set_vertex_z ( pat_mu.vertex().z());
->
+         mu.set_dB ( pat_mu.dB());
->
+         //mu.particleIso ( pat_mu.particleIso());
->
+         mu.set_neutralHadronIso ( pat_mu.neutralHadronIso());
->
+         mu.set_chargedHadronIso ( pat_mu.chargedHadronIso());
->
+         mu.set_trackIso ( pat_mu.trackIso());
->
+         mu.set_photonIso ( pat_mu.photonIso());
->
+         mu.set_puChargedHadronIso ( pat_mu.puChargedHadronIso());
->
+         mu.set_isGlobalMuon ( pat_mu.isGlobalMuon());
->
+         mu.set_isPFMuon ( pat_mu.isPFMuon());
->
+         mu.set_isStandAloneMuon ( pat_mu.isStandAloneMuon());
->
+         mu.set_isTrackerMuon ( pat_mu.isTrackerMuon());
->
+         mu.set_numberOfMatchedStations ( pat_mu.numberOfMatchedStations());
+         reco::TrackRef globalTrack = pat_mu.globalTrack();
+         if(!globalTrack.isNull()){
-<
+           mu.globalTrack_chi2 = globalTrack->chi2();
-<
+           mu.globalTrack_ndof = globalTrack->ndof();
-<
+           mu.globalTrack_d0 = globalTrack->d0();
-<
+           mu.globalTrack_d0Error = globalTrack->d0Error();
-<
+           mu.globalTrack_numberOfValidHits = globalTrack->numberOfValidHits();
-<
+           mu.globalTrack_numberOfLostHits = globalTrack->numberOfLostHits();
->
+           mu.set_globalTrack_chi2 ( globalTrack->chi2());
->
+           mu.set_globalTrack_ndof ( globalTrack->ndof());
->
+           mu.set_globalTrack_d0 ( globalTrack->d0());
->
+           mu.set_globalTrack_d0Error ( globalTrack->d0Error());
->
+           mu.set_globalTrack_numberOfValidHits ( globalTrack->numberOfValidHits());
->
+           mu.set_globalTrack_numberOfLostHits ( globalTrack->numberOfLostHits());
->
+           mu.set_globalTrack_numberOfValidMuonHits(globalTrack->hitPattern().numberOfValidMuonHits() );
+         else{
-<
+           mu.globalTrack_chi2 = 0;
-<
+           mu.globalTrack_ndof = 0;
-<
+           mu.globalTrack_d0 = 0;
-<
+           mu.globalTrack_d0Error = 0;
-<
+           mu.globalTrack_numberOfValidHits = 0;
-<
+           mu.globalTrack_numberOfLostHits = 0;
->
+           mu.set_globalTrack_chi2 ( 0);
->
+           mu.set_globalTrack_ndof ( 0);
->
+           mu.set_globalTrack_d0 ( 0);
->
+           mu.set_globalTrack_d0Error ( 0);
->
+           mu.set_globalTrack_numberOfValidHits ( 0);
->
+           mu.set_globalTrack_numberOfLostHits ( 0);
+         reco::TrackRef innerTrack = pat_mu.innerTrack();
+         if(!innerTrack.isNull()){
-<
+           mu.innerTrack_chi2 = innerTrack->chi2();
-<
+           mu.innerTrack_ndof = innerTrack->ndof();
-<
+           mu.innerTrack_d0 = innerTrack->d0();
-<
+           mu.innerTrack_d0Error = innerTrack->d0Error();
-<
+           mu.innerTrack_numberOfValidHits = innerTrack->numberOfValidHits();
-<
+           mu.innerTrack_numberOfLostHits = innerTrack->numberOfLostHits();
->
+           mu.set_innerTrack_chi2 ( innerTrack->chi2());
->
+           mu.set_innerTrack_ndof ( innerTrack->ndof());
->
+           mu.set_innerTrack_d0 ( innerTrack->d0());
->
+           mu.set_innerTrack_d0Error ( innerTrack->d0Error());
->
+           mu.set_innerTrack_numberOfValidHits ( innerTrack->numberOfValidHits());
->
+           mu.set_innerTrack_numberOfLostHits ( innerTrack->numberOfLostHits());
->
+           mu.set_innerTrack_trackerLayersWithMeasurement ( innerTrack->hitPattern().trackerLayersWithMeasurement());
->
+           mu.set_innerTrack_numberOfValidPixelHits ( innerTrack->hitPattern().numberOfValidPixelHits());
+         else{
-<
+           mu.innerTrack_chi2 = 0;
-<
+           mu.innerTrack_ndof = 0;
-<
+           mu.innerTrack_d0 = 0;
-<
+           mu.innerTrack_d0Error = 0;
-<
+           mu.innerTrack_numberOfValidHits = 0;
-<
+           mu.innerTrack_numberOfLostHits = 0;
->
+           mu.set_innerTrack_chi2 ( 0);
->
+           mu.set_innerTrack_ndof ( 0);
->
+           mu.set_innerTrack_d0 ( 0);
->
+           mu.set_innerTrack_d0Error ( 0);
->
+           mu.set_innerTrack_numberOfValidHits ( 0);
->
+           mu.set_innerTrack_numberOfLostHits ( 0);
->
+           mu.set_innerTrack_trackerLayersWithMeasurement ( 0);
->
+           mu.set_innerTrack_numberOfValidPixelHits ( 0);
+         reco::TrackRef outerTrack = pat_mu.outerTrack();
+         if(!outerTrack.isNull()){
-<
+           mu.outerTrack_chi2 = outerTrack->chi2();
-<
+           mu.outerTrack_ndof = outerTrack->ndof();
-<
+           mu.outerTrack_d0 = outerTrack->d0();
-<
+           mu.outerTrack_d0Error = outerTrack->d0Error();
-<
+           mu.outerTrack_numberOfValidHits = outerTrack->numberOfValidHits();
-<
+           mu.outerTrack_numberOfLostHits = outerTrack->numberOfLostHits();
->
+           mu.set_outerTrack_chi2 ( outerTrack->chi2());
->
+           mu.set_outerTrack_ndof ( outerTrack->ndof());
->
+           mu.set_outerTrack_d0 ( outerTrack->d0());
->
+           mu.set_outerTrack_d0Error ( outerTrack->d0Error());
->
+           mu.set_outerTrack_numberOfValidHits ( outerTrack->numberOfValidHits());
->
+           mu.set_outerTrack_numberOfLostHits ( outerTrack->numberOfLostHits());
+         else{
-<
+           mu.outerTrack_chi2 = 0;
-<
+           mu.outerTrack_ndof = 0;
-<
+           mu.outerTrack_d0 = 0;
-<
+           mu.outerTrack_d0Error = 0;
-<
+           mu.outerTrack_numberOfValidHits = 0;
-<
+           mu.outerTrack_numberOfLostHits = 0;
->
+           mu.set_outerTrack_chi2 ( 0);
->
+           mu.set_outerTrack_ndof ( 0);
->
+           mu.set_outerTrack_d0 ( 0);
->
+           mu.set_outerTrack_d0Error ( 0);
->
+           mu.set_outerTrack_numberOfValidHits ( 0);
->
+           mu.set_outerTrack_numberOfLostHits ( 0);
+         mus[j].push_back(mu);
+         if(fabs(pat_tau.eta()) > tau_etamax) continue;
+         Tau tau;
-<
+         tau.charge =  pat_tau.charge();
-<
+         tau.pt =  pat_tau.pt();
-<
+         tau.eta =  pat_tau.eta();
-<
+         tau.phi =  pat_tau.phi();
-<
+         tau.energy =  pat_tau.energy();
->
+         tau.set_charge( pat_tau.charge());
->
+         tau.set_pt( pat_tau.pt());
->
+         tau.set_eta( pat_tau.eta());
->
+         tau.set_phi( pat_tau.phi());
->
+         tau.set_energy( pat_tau.energy());
->
+         tau.set_decayModeFinding ( pat_tau.tauID("decayModeFinding")>0.5);
->
+         //tau.set_byVLooseCombinedIsolationDeltaBetaCorr  ( pat_tau.tauID("byVLooseCombinedIsolationDeltaBetaCorr")>0.5);
->
+         tau.set_byLooseCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr")>0.5);
->
+         tau.set_byMediumCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr")>0.5);
->
+         tau.set_byTightCombinedIsolationDeltaBetaCorr ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr")>0.5);
->
+         tau.set_byLooseIsolationMVA( pat_tau.tauID("byLooseIsolationMVA")>0.5);
->
+         tau.set_byMediumIsolationMVA( pat_tau.tauID("byMediumIsolationMVA")>0.5);
->
+         tau.set_byTightIsolationMVA( pat_tau.tauID("byTightIsolationMVA")>0.5);
->
+         tau.set_byLooseIsolationMVA2( pat_tau.tauID("byLooseIsolationMVA2")>0.5);
->
+         tau.set_byMediumIsolationMVA2( pat_tau.tauID("byMediumIsolationMVA2")>0.5);
->
+         tau.set_byTightIsolationMVA2( pat_tau.tauID("byTightIsolationMVA2")>0.5);
->
+         tau.set_byLooseCombinedIsolationDeltaBetaCorr3Hits(  pat_tau.tauID("byLooseCombinedIsolationDeltaBetaCorr3Hits")>0.5);
->
+         tau.set_byMediumCombinedIsolationDeltaBetaCorr3Hits ( pat_tau.tauID("byMediumCombinedIsolationDeltaBetaCorr3Hits")>0.5);
->
+         tau.set_byTightCombinedIsolationDeltaBetaCorr3Hits ( pat_tau.tauID("byTightCombinedIsolationDeltaBetaCorr3Hits")>0.5);
->
+         tau.set_againstElectronLooseMVA3  ( pat_tau.tauID("againstElectronLooseMVA3")>0.5);
->
+         tau.set_againstElectronMediumMVA3 ( pat_tau.tauID("againstElectronMediumMVA3")>0.5);
->
+         tau.set_againstElectronTightMVA3 ( pat_tau.tauID("againstElectronTightMVA3")>0.5);
->
+         tau.set_againstElectronVTightMVA3 ( pat_tau.tauID("againstElectronVTightMVA3")>0.5);
->
+         tau.set_againstMuonLoose2 ( pat_tau.tauID("againstMuonLoose2")>0.5);
->
+         tau.set_againstMuonMedium2 ( pat_tau.tauID("againstMuonMedium2")>0.5);
->
+         tau.set_againstMuonTight2 ( pat_tau.tauID("againstMuonTight2")>0.5);
->
+         tau.set_byIsolationMVAraw(  pat_tau.tauID("byIsolationMVAraw"));
->
+         tau.set_byIsolationMVA2raw(  pat_tau.tauID("byIsolationMVA2raw"));
->
+         tau.set_decayMode( pat_tau.decayMode() );
->
+         tau.set_byCombinedIsolationDeltaBetaCorrRaw( pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw"));
->
+         tau.set_byCombinedIsolationDeltaBetaCorrRaw3Hits( pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw3Hits"));
-+
+//       std::cout << pat_tau.tauID("byCombinedIsolationDeltaBetaCorrRaw3Hits") << std::endl;
-+
-+
+//       reco::PFCandidateRef leadPFCand = pat_tau.leadPFCand();
-+
+//       if(!leadPFCand.isNull()){
-+
+//         tau.set_leadPFCand_px ( leadPFCand->px());
-+
+//         tau.set_leadPFCand_py ( leadPFCand->py());
-+
+//         tau.set_leadPFCand_pz ( leadPFCand->pz());
-+
+//       }
-+
+//       else{
-+
+//         tau.set_leadPFCand_px ( 0);
-+
+//         tau.set_leadPFCand_py ( 0);
-+
+//         tau.set_leadPFCand_pz ( 0);
-+
+//       }
+         taus[j].push_back(tau);
-+
+   //-------------- gen jets -------------
-+
-+
+   if(doGenJets){
-+
+     for(size_t j=0; j< genjet_sources.size(); ++j){
-+
-+
+       genjets[j].clear();
-+
-+
+       edm::Handle< std::vector<reco::GenJet> > genjet_handle;
-+
+       iEvent.getByLabel(genjet_sources[j], genjet_handle);
-+
+       const std::vector<reco::GenJet>& gen_jets = *(genjet_handle.product());
-+
-+
+       for (unsigned int i = 0; i < gen_jets.size(); ++i) {
-+
+         const reco::GenJet* gen_jet = &gen_jets[i];
-+
+         if(gen_jet->pt() < genjet_ptmin) continue;
-+
+         if(fabs(gen_jet->eta()) > genjet_etamax) continue;
-+
-+
+         Particle jet;
-+
+         jet.set_charge(gen_jet->charge());
-+
+         jet.set_pt(gen_jet->pt());
-+
+         jet.set_eta(gen_jet->eta());
-+
+         jet.set_phi(gen_jet->phi());
-+
+         jet.set_energy(gen_jet->energy());
-+
-+
+         genjets[j].push_back(jet);
-+
-+
+       }
-+
+     }
-+
+   }
-+
+   // ------------- jets -------------
+   if(doJets){
+     for(size_t j=0; j< jet_sources.size(); ++j){
+//       }
+         Jet jet;
-<
+         jet.charge = pat_jet.charge();
-<
+         jet.pt = pat_jet.pt();
-<
+         jet.eta = pat_jet.eta();
-<
+         jet.phi = pat_jet.phi();
-<
+         jet.energy = pat_jet.energy();
-<
+         jet.numberOfDaughters =pat_jet.numberOfDaughters();
->
+         jet.set_charge(pat_jet.charge());
->
+         jet.set_pt(pat_jet.pt());
->
+         jet.set_eta(pat_jet.eta());
->
+         jet.set_phi(pat_jet.phi());
->
+         jet.set_energy(pat_jet.energy());
->
+         jet.set_numberOfDaughters (pat_jet.numberOfDaughters());
+         const reco::TrackRefVector&  jettracks = pat_jet.associatedTracks();
-<
+         jet.nTracks = jettracks.size();
-<
+         jet.jetArea = pat_jet.jetArea();
-<
+         jet.pileup = pat_jet.pileup();
->
+         jet.set_nTracks ( jettracks.size());
->
+         jet.set_jetArea(pat_jet.jetArea());
->
+         jet.set_pileup(pat_jet.pileup());
+         if(pat_jet.isPFJet()){
-<
+           jet.neutralEmEnergyFraction =pat_jet.neutralEmEnergyFraction();
-<
+           jet.neutralHadronEnergyFraction =pat_jet.neutralHadronEnergyFraction();
-<
+           jet.chargedEmEnergyFraction =pat_jet.chargedEmEnergyFraction();
-<
+           jet.chargedHadronEnergyFraction =pat_jet.chargedHadronEnergyFraction();
-<
+           jet.muonEnergyFraction =pat_jet.muonEnergyFraction();
-<
+           jet.photonEnergyFraction =pat_jet.photonEnergyFraction();
-<
+           jet.chargedMultiplicity =pat_jet.chargedMultiplicity();
-<
+           jet.neutralMultiplicity =pat_jet.neutralMultiplicity();
-<
+           jet.muonMultiplicity =pat_jet.muonMultiplicity();
-<
+           jet.electronMultiplicity =pat_jet.electronMultiplicity();
-<
+           jet.photonMultiplicity =pat_jet.photonMultiplicity();
-<
+         }
-<
-<
+         jecUnc->setJetEta(pat_jet.eta());
-<
+         jecUnc->setJetPt(pat_jet.pt());
-<
+         jet.JEC_uncertainty = jecUnc->getUncertainty(true);
-<
-<
+         jet.btag_simpleSecondaryVertexHighEff=pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
-<
+         jet.btag_simpleSecondaryVertexHighPur=pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
-<
+         jet.btag_combinedSecondaryVertex=pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags");
-<
+         jet.btag_combinedSecondaryVertexMVA=pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
-<
+         jet.btag_jetBProbability=pat_jet.bDiscriminator("jetBProbabilityBJetTags");
-<
+         jet.btag_jetProbability=pat_jet.bDiscriminator("jetProbabilityBJetTags");
->
+           jet.set_neutralEmEnergyFraction (pat_jet.neutralEmEnergyFraction());
->
+           jet.set_neutralHadronEnergyFraction (pat_jet.neutralHadronEnergyFraction());
->
+           jet.set_chargedEmEnergyFraction (pat_jet.chargedEmEnergyFraction());
->
+           jet.set_chargedHadronEnergyFraction (pat_jet.chargedHadronEnergyFraction());
->
+           jet.set_muonEnergyFraction (pat_jet.muonEnergyFraction());
->
+           jet.set_photonEnergyFraction (pat_jet.photonEnergyFraction());
->
+           jet.set_chargedMultiplicity (pat_jet.chargedMultiplicity());
->
+           jet.set_neutralMultiplicity (pat_jet.neutralMultiplicity());
->
+           jet.set_muonMultiplicity (pat_jet.muonMultiplicity());
->
+           jet.set_electronMultiplicity (pat_jet.electronMultiplicity());
->
+           jet.set_photonMultiplicity (pat_jet.photonMultiplicity());
->
+         }
->
+         if(doJECUncertainty){
->
+           jecUnc->setJetEta(pat_jet.eta());
->
+           jecUnc->setJetPt(pat_jet.pt());
->
+           jet.set_JEC_uncertainty(jecUnc->getUncertainty(true));
->
+         }
->
+         jet.set_JEC_factor_raw(pat_jet.jecFactor("Uncorrected"));
->
->
+         jet.set_btag_simpleSecondaryVertexHighEff(pat_jet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
->
+         jet.set_btag_simpleSecondaryVertexHighPur(pat_jet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
->
+         jet.set_btag_combinedSecondaryVertex(pat_jet.bDiscriminator("combinedSecondaryVertexBJetTags"));
->
+         jet.set_btag_combinedSecondaryVertexMVA(pat_jet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
->
+         jet.set_btag_jetBProbability(pat_jet.bDiscriminator("jetBProbabilityBJetTags"));
->
+         jet.set_btag_jetProbability(pat_jet.bDiscriminator("jetProbabilityBJetTags"));
->
->
->
+         const reco::GenJet *genj = pat_jet.genJet();
->
+         if(genj){
-+
+           for(unsigned int k=0; k<genjets->size(); ++k){
-+
+             if(genj->pt()==genjets->at(k).pt() && genj->eta()==genjets->at(k).eta())
-+
+               jet.set_genjet_index(k);
-+
+           }
-+
+//         if( jet.genjet_index()<0){
-+
+//           std::cout<< "genjet not found for " << genj->pt() << "  " << genj->eta() << std::endl;
-+
+//         }
-+
-+
+           if(doAllGenParticles){
-+
+             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
-+
+             for(unsigned int l = 0; l<jetgenps.size(); ++l){
-+
+               for(unsigned int k=0; k< genps.size(); ++k){
-+
+                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
-+
+                   jet.add_genparticles_index(genps[k].index());
-+
+                   break;
-+
+                 }
-+
+               }
-+
+             }
-+
+             if(jet.genparticles_indices().size()!= jetgenps.size())
-+
+               std::cout << "WARNING: Found only " << jet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this jet"<<std::endl;
-+
+           }
-+
-+
+         }
-+
+         jets[j].push_back(jet);
-+
-+
+   // ------------- top jets -------------
+   if(doTopJets){
+     for(size_t j=0; j< topjet_sources.size(); ++j){
+       iEvent.getByLabel(topjet_sources[j], pat_topjets);
+       for (unsigned int i = 0; i < pat_topjets->size(); i++) {
-+
+         const pat::Jet  pat_topjet =  * dynamic_cast<pat::Jet const *>(&pat_topjets->at(i));
+         if(pat_topjet.pt() < topjet_ptmin) continue;
+         if(fabs(pat_topjet.eta()) > topjet_etamax) continue;
+         TopJet topjet;
-<
+         topjet.charge = pat_topjet.charge();
-<
+         topjet.pt = pat_topjet.pt();
-<
+         topjet.eta = pat_topjet.eta();
-<
+         topjet.phi = pat_topjet.phi();
-<
+         topjet.energy = pat_topjet.energy();
-<
+         topjet.numberOfDaughters =pat_topjet.numberOfDaughters();
->
+         topjet.set_charge(pat_topjet.charge());
->
+         topjet.set_pt(pat_topjet.pt());
->
+         topjet.set_eta(pat_topjet.eta());
->
+         topjet.set_phi(pat_topjet.phi());
->
+         topjet.set_energy(pat_topjet.energy());
->
+         topjet.set_numberOfDaughters(pat_topjet.numberOfDaughters());
+         const reco::TrackRefVector&  topjettracks = pat_topjet.associatedTracks();
-<
+         topjet.nTracks = topjettracks.size();
-<
+         topjet.jetArea = pat_topjet.jetArea();
-<
+         topjet.pileup = pat_topjet.pileup();
-<
+//       topjet.neutralEmEnergyFraction =pat_topjet.neutralEmEnergyFraction();
-<
+//       topjet.neutralHadronEnergyFraction =pat_topjet.neutralHadronEnergyFraction();
-<
+//       topjet.chargedEmEnergyFraction =pat_topjet.chargedEmEnergyFraction();
-<
+//       topjet.chargedHadronEnergyFraction =pat_topjet.chargedHadronEnergyFraction();
-<
+//       topjet.muonEnergyFraction =pat_topjet.muonEnergyFraction();
-<
+//       topjet.photonEnergyFraction =pat_topjet.photonEnergyFraction();
-<
+//       topjet.chargedMultiplicity =pat_topjet.chargedMultiplicity();
-<
+//       topjet.neutralMultiplicity =pat_topjet.neutralMultiplicity();
-<
+//       topjet.muonMultiplicity =pat_topjet.muonMultiplicity();
-<
+//       topjet.electronMultiplicity =pat_topjet.electronMultiplicity();
-<
+//       topjet.photonMultiplicity =pat_topjet.photonMultiplicity();
-<
-<
+         jecUnc->setJetEta(pat_topjet.eta());
-<
+         jecUnc->setJetPt(pat_topjet.pt());
-<
+         topjet.JEC_uncertainty = jecUnc->getUncertainty(true);
-<
-<
+         topjet.btag_simpleSecondaryVertexHighEff=pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
-<
+         topjet.btag_simpleSecondaryVertexHighPur=pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
-<
+         topjet.btag_combinedSecondaryVertex=pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags");
-<
+         topjet.btag_combinedSecondaryVertexMVA=pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags");
-<
+         topjet.btag_jetBProbability=pat_topjet.bDiscriminator("jetBProbabilityBJetTags");
-<
+         topjet.btag_jetProbability=pat_topjet.bDiscriminator("jetProbabilityBJetTags");
->
+         topjet.set_nTracks( topjettracks.size());
->
+         topjet.set_jetArea( pat_topjet.jetArea());
->
+         topjet.set_pileup( pat_topjet.pileup());
->
+//       topjet.set_neutralEmEnergyFraction(pat_topjet.neutralEmEnergyFraction());
->
+//       topjet.set_neutralHadronEnergyFraction(pat_topjet.neutralHadronEnergyFraction());
->
+//       topjet.set_chargedEmEnergyFraction(pat_topjet.chargedEmEnergyFraction());
->
+//       topjet.set_chargedHadronEnergyFraction(pat_topjet.chargedHadronEnergyFraction());
->
+//       topjet.set_muonEnergyFraction(pat_topjet.muonEnergyFraction());
->
+//       topjet.set_photonEnergyFraction(pat_topjet.photonEnergyFraction());
->
+//       topjet.set_chargedMultiplicity(pat_topjet.chargedMultiplicity());
->
+//       topjet.set_neutralMultiplicity(pat_topjet.neutralMultiplicity());
->
+//       topjet.set_muonMultiplicity(pat_topjet.muonMultiplicity());
->
+//       topjet.set_electronMultiplicity(pat_topjet.electronMultiplicity());
->
+//       topjet.set_photonMultiplicity(pat_topjet.photonMultiplicity());
->
+         if(doJECUncertainty){
->
+           jecUnc->setJetEta(pat_topjet.eta());
->
+           jecUnc->setJetPt(pat_topjet.pt());
->
+           topjet.set_JEC_uncertainty( jecUnc->getUncertainty(true));
->
+         }
->
+         topjet.set_JEC_factor_raw( pat_topjet.jecFactor("Uncorrected"));
->
->
+         topjet.set_btag_simpleSecondaryVertexHighEff(pat_topjet.bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
->
+         topjet.set_btag_simpleSecondaryVertexHighPur(pat_topjet.bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
->
+         topjet.set_btag_combinedSecondaryVertex(pat_topjet.bDiscriminator("combinedSecondaryVertexBJetTags"));
->
+         topjet.set_btag_combinedSecondaryVertexMVA(pat_topjet.bDiscriminator("combinedSecondaryVertexMVABJetTags"));
->
+         topjet.set_btag_jetBProbability(pat_topjet.bDiscriminator("jetBProbabilityBJetTags"));
->
+         topjet.set_btag_jetProbability(pat_topjet.bDiscriminator("jetProbabilityBJetTags"));
->
->
+         /*
->
+         const reco::GenJet *genj = pat_topjet.genJet();
->
+         if(genj){
->
+           topjet.set_genjet_pt ( genj->pt());
->
+           topjet.set_genjet_eta ( genj->eta());
->
+           topjet.set_genjet_phi ( genj->phi());
->
+           topjet.set_genjet_energy ( genj->energy());
->
+           if(doAllGenParticles){
->
+             std::vector<const reco::GenParticle * > jetgenps = genj->getGenConstituents();
->
+             for(unsigned int l = 0; l<jetgenps.size(); ++l){
->
+               for(unsigned int k=0; k< genps.size(); ++k){
->
+                 if(jetgenps[l]->pt() == genps[k].pt() && jetgenps[l]->pdgId() == genps[k].pdgId()){
->
+                   topjet.add_genparticles_index(genps[k].index());
->
+                   break;
->
+                 }
->
+               }
->
+             }
->
+             if(topjet.genparticles_indices().size()!= jetgenps.size())
->
+               std::cout << "WARNING: Found only " << topjet.genparticles_indices().size() << " from " << jetgenps.size() << " gen particles of this topjet"<<std::endl;
->
+           }
->
+         }
->
+         */
+         for (unsigned int k = 0; k < pat_topjet.numberOfDaughters(); k++) {
+           Particle subjet_v4;
-<
+           subjet_v4.pt = pat_topjet.daughter(k)->p4().pt();
-<
+           subjet_v4.eta = pat_topjet.daughter(k)->p4().eta();
-<
+           subjet_v4.phi = pat_topjet.daughter(k)->p4().phi();
-<
+           subjet_v4.energy = pat_topjet.daughter(k)->p4().E();
-<
+           topjet.subjets.push_back(subjet_v4);
->
->
+           reco::Candidate const * subjetd =  pat_topjet.daughter(k);
->
+           pat::Jet const * patsubjetd = dynamic_cast<pat::Jet const *>(subjetd);
->
+           if(patsubjetd)
->
+           {
->
+              subjet_v4.set_pt(patsubjetd->correctedP4(0).pt());
->
+              subjet_v4.set_eta(patsubjetd->correctedP4(0).eta());
->
+              subjet_v4.set_phi(patsubjetd->correctedP4(0).phi());
->
+              subjet_v4.set_energy(patsubjetd->correctedP4(0).E());
->
+              topjet.add_subjet(subjet_v4);
->
+              //btag info
->
+              topjet.add_subFlavour(patsubjetd->partonFlavour());
->
+              topjet.add_subCSV(patsubjetd->bDiscriminator("combinedSecondaryVertexBJetTags"));
->
+              if (doTagInfos)
->
+                {
->
+                  //ip tag info
->
+                  reco::TaggingVariableList tvlIP=patsubjetd->tagInfoTrackIP("impactParameter")->taggingVariables();
->
+                  topjet.add_subTrackMomentum(tvlIP.getList(reco::btau::trackMomentum,false));
->
+                  topjet.add_subTrackEta(tvlIP.getList(reco::btau::trackEta,false));
->
+                  topjet.add_subTrackEtaRel(tvlIP.getList(reco::btau::trackEtaRel,false));
->
+                  topjet.add_subTrackDeltaR(tvlIP.getList(reco::btau::trackDeltaR,false));
->
+                  topjet.add_subTrackSip3dVal(tvlIP.getList(reco::btau::trackSip3dVal,false));
->
+                  topjet.add_subTrackSip3dSig(tvlIP.getList(reco::btau::trackSip3dSig,false));
->
+                  topjet.add_subTrackSip2dVal(tvlIP.getList(reco::btau::trackSip2dVal,false));
->
+                  topjet.add_subTrackSip2dSig(tvlIP.getList(reco::btau::trackSip2dSig,false));
->
+                  topjet.add_subTrackDecayLenVal(tvlIP.getList(reco::btau::trackDecayLenVal,false));
->
+                  topjet.add_subTrackChi2(tvlIP.getList(reco::btau::trackChi2,false));
->
+                  topjet.add_subTrackNTotalHits(tvlIP.getList(reco::btau::trackNTotalHits,false));
->
+                  topjet.add_subTrackNPixelHits(tvlIP.getList(reco::btau::trackNPixelHits,false));
->
+                  topjet.add_subTrackPtRel(tvlIP.getList(reco::btau::trackPtRel,false));
->
+                  topjet.add_subTrackPPar(tvlIP.getList(reco::btau::trackPPar,false));
->
+                  topjet.add_subTrackPtRatio(tvlIP.getList(reco::btau::trackPtRatio,false));
->
+                  topjet.add_subTrackPParRatio(tvlIP.getList(reco::btau::trackPParRatio,false));
->
+                  topjet.add_subTrackJetDistVal(tvlIP.getList(reco::btau::trackJetDistVal,false));
->
+                  topjet.add_subTrackJetDistSig(tvlIP.getList(reco::btau::trackJetDistSig,false));
->
+                  topjet.add_subTrackGhostTrackDistVal(tvlIP.getList(reco::btau::trackGhostTrackDistVal,false));
->
+                  topjet.add_subTrackGhostTrackDistSig(tvlIP.getList(reco::btau::trackGhostTrackDistSig,false));
->
+                  topjet.add_subTrackGhostTrackWeight(tvlIP.getList(reco::btau::trackGhostTrackWeight,false));
->
+                  //sv tag info
->
+                  reco::TaggingVariableList tvlSV=patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->taggingVariables();
->
+                  topjet.add_subFlightDistance2dVal(tvlSV.getList(reco::btau::flightDistance2dVal,false));
->
+                  topjet.add_subFlightDistance2dSig(tvlSV.getList(reco::btau::flightDistance2dSig,false));
->
+                  topjet.add_subFlightDistance3dVal(tvlSV.getList(reco::btau::flightDistance3dVal,false));
->
+                  topjet.add_subFlightDistance3dSig(tvlSV.getList(reco::btau::flightDistance3dSig,false));
->
+                  topjet.add_subVertexJetDeltaR(tvlSV.getList(reco::btau::vertexJetDeltaR,false));
->
+                  topjet.add_subJetNSecondaryVertices(tvlSV.get(reco::btau::jetNSecondaryVertices,-9999));
->
+                  topjet.add_subVertexNTracks(tvlSV.get(reco::btau::vertexNTracks,-9999));
->
+                  std::vector<TLorentzVector> vp4; vp4.clear();
->
+                  std::vector<float> vchi2; vchi2.clear();
->
+                  std::vector<float> vndof; vndof.clear();
->
+                  std::vector<float> vchi2ndof; vchi2ndof.clear();
->
+                  std::vector<float> vtrsize; vtrsize.clear();
->
+                  for(unsigned int i=0; i<patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->nVertices(); i++)
->
+                    {
->
+                      ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > p4 = patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).p4();
->
+                      vp4.push_back(TLorentzVector(p4.px(),p4.py(),p4.pz(),p4.e()));
->
+                      vchi2.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).chi2());
->
+                      vndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).ndof());
->
+                      vchi2ndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).normalizedChi2());
->
+                      vtrsize.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).tracksSize());
->
+                    }
->
+                  topjet.add_subSecondaryVertex(vp4);
->
+                  topjet.add_subVertexChi2(vchi2);
->
+                  topjet.add_subVertexNdof(vndof);
->
+                  topjet.add_subVertexNormalizedChi2(vchi2ndof);
->
+                  topjet.add_subVertexTracksSize(vtrsize);
->
+                  //try computer
->
+                  edm::ESHandle<JetTagComputer> computerHandle;
->
+                  iSetup.get<JetTagComputerRecord>().get( SVComputer_.label(), computerHandle );
->
+                  const GenericMVAJetTagComputer *computer = dynamic_cast<const GenericMVAJetTagComputer*>( computerHandle.product() );
->
+                  if(computer)
->
+                    {
->
+                      computer->passEventSetup(iSetup);
->
+                      std::vector<const reco::BaseTagInfo*>  baseTagInfos;
->
+                      baseTagInfos.push_back(patsubjetd->tagInfoTrackIP("impactParameter") );
->
+                      baseTagInfos.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex") );
->
+                      JetTagComputer::TagInfoHelper helper(baseTagInfos);
->
+                      reco::TaggingVariableList vars = computer->taggingVariables(helper);
->
+                      topjet.add_subVertexMassJTC(vars.get(reco::btau::vertexMass,-9999));
->
+                      topjet.add_subVertexCategoryJTC(vars.get(reco::btau::vertexCategory,-9999));
->
+                      topjet.add_subVertexEnergyRatioJTC(vars.get(reco::btau::vertexEnergyRatio,-9999));
->
+                      topjet.add_subTrackSip3dSigAboveCharmJTC(vars.get(reco::btau::trackSip3dSigAboveCharm,-9999));
->
+                    }
->
+                }
->
+           }
->
+           else
->
+             {
->
+               //filling only standard information in case the subjet has not been pat-tified during the pattuples production
->
+               subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt());
->
+               subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta());
->
+               subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi());
->
+               subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E());
->
+               topjet.add_subjet(subjet_v4);
->
+             }
->
->
+         topjets[j].push_back(topjet);
-+
-+
-+
+   // ------------- generator top jets -------------
-+
+   if(doGenTopJets){
-+
+     for(size_t j=0; j< gentopjet_sources.size(); ++j){
-+
-+
+       gentopjets[j].clear();
-+
-+
+       edm::Handle<reco::BasicJetCollection> reco_gentopjets;
-+
+       //edm::Handle<std::vector<reco::Jet> > reco_gentopjets;
-+
+       iEvent.getByLabel(gentopjet_sources[j], reco_gentopjets);
-+
-+
+       for (unsigned int i = 0; i < reco_gentopjets->size(); i++) {
-+
-+
+         const reco::BasicJet  reco_gentopjet =  reco_gentopjets->at(i);
-+
+         if(reco_gentopjet.pt() < gentopjet_ptmin) continue;
-+
+         if(fabs(reco_gentopjet.eta()) > gentopjet_etamax) continue;
-+
-+
+         TopJet gentopjet;
-+
+         gentopjet.set_charge(reco_gentopjet.charge());
-+
+         gentopjet.set_pt(reco_gentopjet.pt());
-+
+         gentopjet.set_eta(reco_gentopjet.eta());
-+
+         gentopjet.set_phi(reco_gentopjet.phi());
-+
+         gentopjet.set_energy(reco_gentopjet.energy());
-+
+         gentopjet.set_numberOfDaughters(reco_gentopjet.numberOfDaughters());
-+
-+
+         for (unsigned int k = 0; k < reco_gentopjet.numberOfDaughters(); k++) {
-+
+           Particle subjet_v4;
-+
+           subjet_v4.set_pt(reco_gentopjet.daughter(k)->p4().pt());
-+
+           subjet_v4.set_eta(reco_gentopjet.daughter(k)->p4().eta());
-+
+           subjet_v4.set_phi(reco_gentopjet.daughter(k)->p4().phi());
-+
+           subjet_v4.set_energy(reco_gentopjet.daughter(k)->p4().E());
-+
+           gentopjet.add_subjet(subjet_v4);
-+
+         }
-+
+         gentopjets[j].push_back(gentopjet);
-+
+       }
-+
+     }
-+
+   }
+   // ------------- photons -------------
+   if(doPhotons){
+       for (unsigned int i = 0; i < pat_photons.size(); ++i) {
+         pat::Photon pat_photon = pat_photons[i];
+         Photon ph;
-<
+         ph.charge = 0;
-<
+         ph.pt =  pat_photon.pt();
-<
+         ph.eta =  pat_photon.eta();
-<
+         ph.phi =  pat_photon.phi();
-<
+         ph.energy =  pat_photon.energy();
-<
+         ph.vertex_x = pat_photon.vertex().x();
-<
+         ph.vertex_y = pat_photon.vertex().y();
-<
+         ph.vertex_z = pat_photon.vertex().z();
-<
+         ph.supercluster_eta = pat_photon.superCluster()->eta();
-<
+         ph.supercluster_phi = pat_photon.superCluster()->phi();
-<
+//       ph.neutralHadronIso = pat_photon.neutralHadronIso();
-<
+//       ph.chargedHadronIso = pat_photon.chargedHadronIso();
-<
+         ph.trackIso = pat_photon.trackIso();
->
+         ph.set_charge(0);
->
+         ph.set_pt( pat_photon.pt());
->
+         ph.set_eta( pat_photon.eta());
->
+         ph.set_phi( pat_photon.phi());
->
+         ph.set_energy( pat_photon.energy());
->
+         ph.set_vertex_x(pat_photon.vertex().x());
->
+         ph.set_vertex_y(pat_photon.vertex().y());
->
+         ph.set_vertex_z(pat_photon.vertex().z());
->
+         ph.set_supercluster_eta(pat_photon.superCluster()->eta());
->
+         ph.set_supercluster_phi(pat_photon.superCluster()->phi());
->
+//       ph.set_neutralHadronIso(pat_photon.neutralHadronIso());
->
+//       ph.set_chargedHadronIso(pat_photon.chargedHadronIso());
->
+         ph.set_trackIso(pat_photon.trackIso());
+         phs[j].push_back(ph);
+       else{
+         pat::MET pat_met = pat_mets[0];
-<
+         met[j].pt=pat_met.pt();
-<
+         met[j].phi=pat_met.phi();
-<
+         met[j].mEtSig= pat_met.mEtSig();
->
+         met[j].set_pt(pat_met.pt());
->
+         met[j].set_phi(pat_met.phi());
->
+         met[j].set_mEtSig(pat_met.mEtSig());
+     triggerResults.clear();
+     triggerNames.clear();
-<
+     L1_prescale.clear();
-<
+     HLT_prescale.clear();
->
+//      L1_prescale.clear();
->
+//      HLT_prescale.clear();
+     edm::Service<edm::service::TriggerNamesService> tns;
+     std::vector<std::string> triggerNames_all;
+       //triggerResults.insert(std::pair<std::string, bool>(triggerNames[i],trig.accept(i)));
+       triggerResults.push_back(trig.accept(i));
+       if(newrun) triggerNames.push_back(triggerNames_all[i]);
-<
+       if(isRealData){
-<
+         std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
-<
+         L1_prescale.push_back(pre.first);
-<
+         HLT_prescale.push_back(pre.second);
-<
+       }
->
+//        if(isRealData){
->
+//       std::pair<int, int> pre=hlt_cfg.prescaleValues(iEvent, iSetup, triggerNames_all[i]);
->
+//       L1_prescale.push_back(pre.first);
->
+//       HLT_prescale.push_back(pre.second);
->
+//       //std::cout <<  triggerNames_all[i] << " " << pre.first << " " <<pre.second << "   " << hlt_cfg.prescaleValue(iEvent, iSetup, triggerNames_all[i]) << std::endl;
->
+//        }
+     //    for(std::map<std::string, bool>::const_iterator iter = triggerResults.begin(); iter!=triggerResults.end(); iter++){
+     //      std::cout << (*iter).first << "   " << (*iter).second << std::endl;
+     newrun=false;
-–
+   // ------------- generator info -------------
-–
-–
+   if(doGenInfo){
-–
+     genInfo.weights.clear();
-–
+     genInfo.binningValues.clear();
-–
+     genps.clear();
-–
-–
+     edm::Handle<GenEventInfoProduct> genEventInfoProduct;
-–
+     iEvent.getByLabel("generator", genEventInfoProduct);
-–
+     const GenEventInfoProduct& genEventInfo = *(genEventInfoProduct.product());
-–
-–
+     genInfo.binningValues = genEventInfo.binningValues();
-–
+     genInfo.weights = genEventInfo.weights();
-–
+     genInfo.alphaQCD = genEventInfo.alphaQCD();
-–
+     genInfo.alphaQED = genEventInfo.alphaQED();
-–
+     genInfo.qScale = genEventInfo.qScale();
-–
-–
+     const gen::PdfInfo* pdf = genEventInfo.pdf();
-–
+     if(pdf){
-–
+       genInfo.pdf_id1=pdf->id.first;
-–
+       genInfo.pdf_id2=pdf->id.second;
-–
+       genInfo.pdf_x1=pdf->x.first;
-–
+       genInfo.pdf_x2=pdf->x.second;
-–
+       genInfo.pdf_scalePDF=pdf->scalePDF;
-–
+       genInfo.pdf_xPDF1=pdf->xPDF.first;
-–
+       genInfo.pdf_xPDF2=pdf->xPDF.second;
-–
+     }
-–
+     else{
-–
+       genInfo.pdf_id1=-999;
-–
+       genInfo.pdf_id2=-999;
-–
+       genInfo.pdf_x1=-999;
-–
+       genInfo.pdf_x2=-999;
-–
+       genInfo.pdf_scalePDF=-999;
-–
+       genInfo.pdf_xPDF1=-999;
-–
+       genInfo.pdf_xPDF2=-999;
-–
+     }
-–
-–
+     edm::Handle<std::vector<PileupSummaryInfo> > pus;
-–
+     iEvent.getByLabel(edm::InputTag("addPileupInfo"), pus);
-–
+     genInfo.pileup_NumInteractions_intime=0;
-–
+     genInfo.pileup_NumInteractions_ootbefore=0;
-–
+     genInfo.pileup_NumInteractions_ootafter=0;
-–
+     if(pus.isValid()){
-–
+       genInfo.pileup_TrueNumInteractions = (float) pus->at(0).getTrueNumInteractions();
-–
+       for(size_t i=0; i<pus->size(); ++i){
-–
+         if(pus->at(i).getBunchCrossing() == 0) // intime pileup
-–
+            genInfo.pileup_NumInteractions_intime += pus->at(i).getPU_NumInteractions();
-–
+         else if(pus->at(i).getBunchCrossing() == -1){ // oot pileup before
-–
+            genInfo.pileup_NumInteractions_ootbefore += pus->at(i).getPU_NumInteractions();
-–
+         }
-–
+         else if(pus->at(i).getBunchCrossing() == +1){ // oot pileup before
-–
+           genInfo.pileup_NumInteractions_ootafter += pus->at(i).getPU_NumInteractions();
-–
+         }
-–
+       }
-–
+     }
-–
-–
+     edm::Handle<reco::GenParticleCollection> genPartColl;
-–
+     iEvent.getByLabel(edm::InputTag("genParticles"), genPartColl);
-–
+     int index=-1;
-–
+     for(reco::GenParticleCollection::const_iterator iter = genPartColl->begin(); iter != genPartColl->end(); ++ iter){
-–
+       index++;
-–
-–
+       //write out only top quarks and status 3 particles (works fine only for MadGraph)
-–
+       if(abs(iter->pdgId())==6 || iter->status()==3){
-–
+         GenParticle genp;
-–
+         genp.charge = iter->charge();;
-–
+         genp.pt = iter->p4().pt();
-–
+         genp.eta = iter->p4().eta();
-–
+         genp.phi = iter->p4().phi();
-–
+         genp.energy = iter->p4().E();
-–
+         genp.index =index;
-–
+         genp.status = iter->status();
-–
+         genp.pdgId = iter->pdgId();
-–
-–
+         genp.mother1=-1;
-–
+         genp.mother2=-1;
-–
+         genp.daughter1=-1;
-–
+         genp.daughter2=-1;
-–
-–
+         int nm=iter->numberOfMothers();
-–
+         int nd=iter->numberOfDaughters();
-–
-–
+         if (nm>0) genp.mother1 = iter->motherRef(0).key();
-–
+         if (nm>1) genp.mother2 = iter->motherRef(nm-1).key();
-–
+         if (nd>0) genp.daughter1 = iter->daughterRef(0).key();
-–
+         if (nd>1) genp.daughter2 = iter->daughterRef(nd-1).key();
-–
-–
+         genps.push_back(genp);
-–
+       }
-–
+     }
-–
-–
+   }
-–
+   tr->Fill();
-<
->
+   if(doLumiInfo)
->
+     previouslumiblockwasfilled=true;
+void
+NtupleWriter::beginJob()
-+
+  if(doLumiInfo){
-+
+    totalRecLumi=0;
-+
+    totalDelLumi=0;
-+
+    previouslumiblockwasfilled=false;
-+
+  }
+// ------------ method called once each job just after ending the event loop  ------------
+void
+NtupleWriter::beginRun(edm::Run const& iRun, edm::EventSetup const&  iSetup)
-<
+  bool setup_changed = false;
-<
+  hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
-<
+  newrun=true;
-<
-<
+  edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
-<
+  iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
-<
+  JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
-<
+  jecUnc = new JetCorrectionUncertainty(JetCorPar);
->
+  if(doTrigger){
->
+    //bool setup_changed = false;
->
+    //hlt_cfg.init(iRun, iSetup, "HLT", setup_changed);
->
+    newrun=true;
->
+  }
-+
+  if(doJets || doTopJets){
-+
+    if(doJECUncertainty){
-+
+      edm::ESHandle<JetCorrectorParametersCollection> JetCorParColl;
-+
+      iSetup.get<JetCorrectionsRecord>().get("AK5PF",JetCorParColl);
-+
+      JetCorrectorParameters const & JetCorPar = (*JetCorParColl)["Uncertainty"];
-+
+      jecUnc = new JetCorrectionUncertainty(JetCorPar);
-+
+    }
-+
+  }
+// ------------ method called when ending the processing of a run  ------------
+void
+NtupleWriter::endRun(edm::Run const&, edm::EventSetup const&)
-+
+  if(doLumiInfo)
-+
+    std::cout << "total integ. luminosity: " << totalDelLumi <<"(del) " << totalRecLumi << "(rec)" << std::endl;
+// ------------ method called when starting to processes a luminosity block  ------------
+void
-<
+NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&)
->
+NtupleWriter::beginLuminosityBlock(edm::LuminosityBlock const& lumi, edm::EventSetup const&)
-+
+  if(doLumiInfo){
-+
+    edm::Handle<LumiSummary> l;
-+
+    lumi.getByLabel("lumiProducer", l);
-+
-+
+    //add lumi of lumi blocks without any event to next lumiblock
-+
+    if(previouslumiblockwasfilled){
-+
+      intgRecLumi=0;
-+
+      intgDelLumi=0;
-+
+    }
-+
+    previouslumiblockwasfilled=false;
-+
-+
+    if (l.isValid()){;
-+
+      intgRecLumi+=l->intgRecLumi()*6.37;
-+
+      intgDelLumi+=l->intgDelLumi()*6.37;
-+
+      totalRecLumi+=l->intgRecLumi()*6.37;
-+
+      totalDelLumi+=l->intgDelLumi()*6.37;
-+
+    }
-+
+    //std::cout << "this lb: " <<l->intgRecLumi()*6.37 <<"   " << l->intgDelLumi()*6.37<<std::endl;
-+
+    //std::cout << "summed:  "<< intgRecLumi << "   " << intgDelLumi << std::endl;
-+
+  }
+// ------------ method called when ending the processing of a luminosity block  ------------

Diff Legend

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

Comparing UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc (file contents): Revision 1.6 by peiffer, Wed Apr 11 15:15:44 2012 UTC vs. Revision 1.26 by peiffer, Mon Jun 10 13:33:05 2013 UTC

Diff Legend

Comparing UserCode/UHHAnalysis/NtupleWriter/src/NtupleWriter.cc (file contents):
Revision 1.6 by peiffer, Wed Apr 11 15:15:44 2012 UTC vs.
Revision 1.26 by peiffer, Mon Jun 10 13:33:05 2013 UTC