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 |
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"); |
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); |
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"); |
641 |
|
} |
642 |
|
} |
643 |
|
|
644 |
+ |
//-------------- gen jets ------------- |
645 |
+ |
|
646 |
+ |
if(doGenJets){ |
647 |
+ |
for(size_t j=0; j< genjet_sources.size(); ++j){ |
648 |
+ |
|
649 |
+ |
genjets[j].clear(); |
650 |
+ |
|
651 |
+ |
edm::Handle< std::vector<reco::GenJet> > genjet_handle; |
652 |
+ |
iEvent.getByLabel(genjet_sources[j], genjet_handle); |
653 |
+ |
const std::vector<reco::GenJet>& gen_jets = *(genjet_handle.product()); |
654 |
+ |
|
655 |
+ |
for (unsigned int i = 0; i < gen_jets.size(); ++i) { |
656 |
+ |
pat::Jet gen_jet = gen_jets[i]; |
657 |
+ |
if(gen_jet.pt() < genjet_ptmin) continue; |
658 |
+ |
if(fabs(gen_jet.eta()) > genjet_etamax) continue; |
659 |
+ |
|
660 |
+ |
Particle jet; |
661 |
+ |
jet.set_charge(gen_jet.charge()); |
662 |
+ |
jet.set_pt(gen_jet.pt()); |
663 |
+ |
jet.set_eta(gen_jet.eta()); |
664 |
+ |
jet.set_phi(gen_jet.phi()); |
665 |
+ |
jet.set_energy(gen_jet.energy()); |
666 |
+ |
|
667 |
+ |
genjets[j].push_back(jet); |
668 |
+ |
|
669 |
+ |
} |
670 |
+ |
} |
671 |
+ |
} |
672 |
+ |
|
673 |
|
// ------------- top jets ------------- |
674 |
|
if(doTopJets){ |
675 |
|
for(size_t j=0; j< topjet_sources.size(); ++j){ |
747 |
|
|
748 |
|
for (unsigned int k = 0; k < pat_topjet.numberOfDaughters(); k++) { |
749 |
|
Particle subjet_v4; |
750 |
< |
subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt()); |
751 |
< |
subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta()); |
752 |
< |
subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi()); |
753 |
< |
subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E()); |
754 |
< |
topjet.add_subjet(subjet_v4); |
750 |
> |
|
751 |
> |
reco::Candidate const * subjetd = pat_topjet.daughter(k); |
752 |
> |
pat::Jet const * patsubjetd = dynamic_cast<pat::Jet const *>(subjetd); |
753 |
> |
if(patsubjetd) |
754 |
> |
{ |
755 |
> |
subjet_v4.set_pt(patsubjetd->correctedP4(0).pt()); |
756 |
> |
subjet_v4.set_eta(patsubjetd->correctedP4(0).eta()); |
757 |
> |
subjet_v4.set_phi(patsubjetd->correctedP4(0).phi()); |
758 |
> |
subjet_v4.set_energy(patsubjetd->correctedP4(0).E()); |
759 |
> |
topjet.add_subjet(subjet_v4); |
760 |
> |
//btag info |
761 |
> |
topjet.add_subFlavour(patsubjetd->partonFlavour()); |
762 |
> |
topjet.add_subCSV(patsubjetd->bDiscriminator("combinedSecondaryVertexBJetTags")); |
763 |
> |
if (doTagInfos) |
764 |
> |
{ |
765 |
> |
//ip tag info |
766 |
> |
reco::TaggingVariableList tvlIP=patsubjetd->tagInfoTrackIP("impactParameter")->taggingVariables(); |
767 |
> |
topjet.add_subTrackMomentum(tvlIP.getList(reco::btau::trackMomentum,false)); |
768 |
> |
topjet.add_subTrackEta(tvlIP.getList(reco::btau::trackEta,false)); |
769 |
> |
topjet.add_subTrackEtaRel(tvlIP.getList(reco::btau::trackEtaRel,false)); |
770 |
> |
topjet.add_subTrackDeltaR(tvlIP.getList(reco::btau::trackDeltaR,false)); |
771 |
> |
topjet.add_subTrackSip3dVal(tvlIP.getList(reco::btau::trackSip3dVal,false)); |
772 |
> |
topjet.add_subTrackSip3dSig(tvlIP.getList(reco::btau::trackSip3dSig,false)); |
773 |
> |
topjet.add_subTrackSip2dVal(tvlIP.getList(reco::btau::trackSip2dVal,false)); |
774 |
> |
topjet.add_subTrackSip2dSig(tvlIP.getList(reco::btau::trackSip2dSig,false)); |
775 |
> |
topjet.add_subTrackDecayLenVal(tvlIP.getList(reco::btau::trackDecayLenVal,false)); |
776 |
> |
topjet.add_subTrackChi2(tvlIP.getList(reco::btau::trackChi2,false)); |
777 |
> |
topjet.add_subTrackNTotalHits(tvlIP.getList(reco::btau::trackNTotalHits,false)); |
778 |
> |
topjet.add_subTrackNPixelHits(tvlIP.getList(reco::btau::trackNPixelHits,false)); |
779 |
> |
topjet.add_subTrackPtRel(tvlIP.getList(reco::btau::trackPtRel,false)); |
780 |
> |
topjet.add_subTrackPPar(tvlIP.getList(reco::btau::trackPPar,false)); |
781 |
> |
topjet.add_subTrackPtRatio(tvlIP.getList(reco::btau::trackPtRatio,false)); |
782 |
> |
topjet.add_subTrackPParRatio(tvlIP.getList(reco::btau::trackPParRatio,false)); |
783 |
> |
topjet.add_subTrackJetDistVal(tvlIP.getList(reco::btau::trackJetDistVal,false)); |
784 |
> |
topjet.add_subTrackJetDistSig(tvlIP.getList(reco::btau::trackJetDistSig,false)); |
785 |
> |
topjet.add_subTrackGhostTrackDistVal(tvlIP.getList(reco::btau::trackGhostTrackDistVal,false)); |
786 |
> |
topjet.add_subTrackGhostTrackDistSig(tvlIP.getList(reco::btau::trackGhostTrackDistSig,false)); |
787 |
> |
topjet.add_subTrackGhostTrackWeight(tvlIP.getList(reco::btau::trackGhostTrackWeight,false)); |
788 |
> |
//sv tag info |
789 |
> |
reco::TaggingVariableList tvlSV=patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->taggingVariables(); |
790 |
> |
topjet.add_subFlightDistance2dVal(tvlSV.getList(reco::btau::flightDistance2dVal,false)); |
791 |
> |
topjet.add_subFlightDistance2dSig(tvlSV.getList(reco::btau::flightDistance2dSig,false)); |
792 |
> |
topjet.add_subFlightDistance3dVal(tvlSV.getList(reco::btau::flightDistance3dVal,false)); |
793 |
> |
topjet.add_subFlightDistance3dSig(tvlSV.getList(reco::btau::flightDistance3dSig,false)); |
794 |
> |
topjet.add_subVertexJetDeltaR(tvlSV.getList(reco::btau::vertexJetDeltaR,false)); |
795 |
> |
topjet.add_subJetNSecondaryVertices(tvlSV.get(reco::btau::jetNSecondaryVertices,-9999)); |
796 |
> |
topjet.add_subVertexNTracks(tvlSV.get(reco::btau::vertexNTracks,-9999)); |
797 |
> |
std::vector<TLorentzVector> vp4; vp4.clear(); |
798 |
> |
std::vector<float> vchi2; vchi2.clear(); |
799 |
> |
std::vector<float> vndof; vndof.clear(); |
800 |
> |
std::vector<float> vchi2ndof; vchi2ndof.clear(); |
801 |
> |
std::vector<float> vtrsize; vtrsize.clear(); |
802 |
> |
for(unsigned int i=0; i<patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->nVertices(); i++) |
803 |
> |
{ |
804 |
> |
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > p4 = patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).p4(); |
805 |
> |
vp4.push_back(TLorentzVector(p4.px(),p4.py(),p4.pz(),p4.e())); |
806 |
> |
vchi2.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).chi2()); |
807 |
> |
vndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).ndof()); |
808 |
> |
vchi2ndof.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).normalizedChi2()); |
809 |
> |
vtrsize.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex")->secondaryVertex(i).tracksSize()); |
810 |
> |
} |
811 |
> |
topjet.add_subSecondaryVertex(vp4); |
812 |
> |
topjet.add_subVertexChi2(vchi2); |
813 |
> |
topjet.add_subVertexNdof(vndof); |
814 |
> |
topjet.add_subVertexNormalizedChi2(vchi2ndof); |
815 |
> |
topjet.add_subVertexTracksSize(vtrsize); |
816 |
> |
//try computer |
817 |
> |
edm::ESHandle<JetTagComputer> computerHandle; |
818 |
> |
iSetup.get<JetTagComputerRecord>().get( SVComputer_.label(), computerHandle ); |
819 |
> |
const GenericMVAJetTagComputer *computer = dynamic_cast<const GenericMVAJetTagComputer*>( computerHandle.product() ); |
820 |
> |
if(computer) |
821 |
> |
{ |
822 |
> |
computer->passEventSetup(iSetup); |
823 |
> |
std::vector<const reco::BaseTagInfo*> baseTagInfos; |
824 |
> |
baseTagInfos.push_back(patsubjetd->tagInfoTrackIP("impactParameter") ); |
825 |
> |
baseTagInfos.push_back(patsubjetd->tagInfoSecondaryVertex("secondaryVertex") ); |
826 |
> |
JetTagComputer::TagInfoHelper helper(baseTagInfos); |
827 |
> |
reco::TaggingVariableList vars = computer->taggingVariables(helper); |
828 |
> |
topjet.add_subVertexMassJTC(vars.get(reco::btau::vertexMass,-9999)); |
829 |
> |
topjet.add_subVertexCategoryJTC(vars.get(reco::btau::vertexCategory,-9999)); |
830 |
> |
topjet.add_subVertexEnergyRatioJTC(vars.get(reco::btau::vertexEnergyRatio,-9999)); |
831 |
> |
topjet.add_subTrackSip3dSigAboveCharmJTC(vars.get(reco::btau::trackSip3dSigAboveCharm,-9999)); |
832 |
> |
} |
833 |
> |
} |
834 |
> |
} |
835 |
> |
else |
836 |
> |
{ |
837 |
> |
//filling only standard information in case the subjet has not been pat-tified during the pattuples production |
838 |
> |
subjet_v4.set_pt(pat_topjet.daughter(k)->p4().pt()); |
839 |
> |
subjet_v4.set_eta(pat_topjet.daughter(k)->p4().eta()); |
840 |
> |
subjet_v4.set_phi(pat_topjet.daughter(k)->p4().phi()); |
841 |
> |
subjet_v4.set_energy(pat_topjet.daughter(k)->p4().E()); |
842 |
> |
topjet.add_subjet(subjet_v4); |
843 |
> |
} |
844 |
> |
|
845 |
> |
|
846 |
|
} |
847 |
|
topjets[j].push_back(topjet); |
848 |
|
} |
849 |
|
} |
850 |
|
} |
851 |
< |
|
852 |
< |
|
851 |
> |
|
852 |
> |
|
853 |
|
// ------------- generator top jets ------------- |
854 |
|
if(doGenTopJets){ |
855 |
|
for(size_t j=0; j< gentopjet_sources.size(); ++j){ |
856 |
|
|
857 |
|
gentopjets[j].clear(); |
858 |
< |
|
858 |
> |
|
859 |
|
edm::Handle<reco::BasicJetCollection> reco_gentopjets; |
860 |
|
//edm::Handle<std::vector<reco::Jet> > reco_gentopjets; |
861 |
|
iEvent.getByLabel(gentopjet_sources[j], reco_gentopjets); |