ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitProd/TreeFiller/src/FillerPFJets.cc
Revision: 1.11
Committed: Fri Mar 26 14:19:00 2010 UTC (15 years, 1 month ago) by sixie
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_018, Mit_017, Mit_017pre3, Mit_017pre2, Mit_017pre1, Mit_016, Mit_015b, Mit_015a, Mit_015, Mit_014e, Mit_014d, Mit_014c, Mit_014b, Mit_014a, Mit_014, Mit_014pre3, Mit_014pre2, Mit_014pre1, Mit_013d, Mit_013c, Mit_013b, Mit_013a
Changes since 1.10: +7 -2 lines
Log Message:
Add jet moments

File Contents

# User Rev Content
1 sixie 1.11 // $Id: FillerPFJets.cc,v 1.10 2010/03/18 20:21:01 bendavid Exp $
2 bendavid 1.1
3     #include "MitProd/TreeFiller/interface/FillerPFJets.h"
4 loizides 1.5 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
5     #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
6     #include "DataFormats/JetReco/interface/PFJet.h"
7     #include "DataFormats/BTauReco/interface/JetTag.h"
8 bendavid 1.1 #include "SimDataFormats/JetMatching/interface/JetFlavour.h"
9     #include "SimDataFormats/JetMatching/interface/JetFlavourMatching.h"
10     #include "SimDataFormats/JetMatching/interface/MatchedPartons.h"
11     #include "SimDataFormats/JetMatching/interface/JetMatchedPartons.h"
12     #include "JetMETCorrections/Objects/interface/JetCorrector.h"
13 loizides 1.5 #include "MitAna/DataTree/interface/Names.h"
14     #include "MitAna/DataTree/interface/PFJetCol.h"
15     #include "MitProd/ObjectService/interface/ObjectService.h"
16 bendavid 1.1
17     using namespace std;
18     using namespace edm;
19     using namespace mithep;
20    
21     //--------------------------------------------------------------------------------------------------
22     FillerPFJets::FillerPFJets(const ParameterSet &cfg, const char *name, bool active) :
23     BaseFiller(cfg,name,active),
24     flavorMatchingActive_(Conf().getUntrackedParameter<bool>("flavorMatchingActive",true)),
25     bTaggingActive_(Conf().getUntrackedParameter<bool>("bTaggingActive",true)),
26     jetToVertexActive_(Conf().getUntrackedParameter<bool>("jetToVertexActive",true)),
27     jetCorrectionsActive_(Conf().getUntrackedParameter<bool>("jetCorrectionsActive",true)),
28     edmName_(Conf().getUntrackedParameter<string>("edmName","recoPFJets:iterativeCone5PFJets")),
29     mitName_(Conf().getUntrackedParameter<string>("mitName","ItrCone5PFJets")),
30     jetToVertexAlphaName_(Conf().getUntrackedParameter<string>
31     ("jetToVertexAlphaName","jetToVertexAlpha")),
32     jetToVertexBetaName_(Conf().getUntrackedParameter<string>
33     ("jetToVertexBetaName","jetToVertexBetaName")),
34     L2JetCorrectorName_(Conf().getUntrackedParameter<string>
35     ("L2JetCorrectorName","L2JetCorrectorName")),
36     L3JetCorrectorName_(Conf().getUntrackedParameter<string>
37     ("L3JetCorrectorName","L3JetCorrectorName")),
38     flavorMatchingByReferenceName_(Conf().getUntrackedParameter<string>
39     ("flavorMatchingByReferenceName","srcByReference")),
40     flavorMatchingDefinition_(Conf().getUntrackedParameter<string>
41     ("flavorMatchingDefinition","Algorithmic")),
42     jetProbabilityBJetTagsName_(Conf().getUntrackedParameter<string>
43     ("JetProbabilityBJetTagsName","jetProbabilityBJetTags")),
44     jetBProbabilityBJetTagsName_(Conf().getUntrackedParameter<string>
45     ("JetBProbabilityBJetTagsName","jetBProbabilityBJetTags")),
46     simpleSecondaryVertexBJetTagsName_(Conf().getUntrackedParameter<string>
47     ("SimpleSecondaryVertexBJetTagsName","simpleSecondaryVertexBJetTags")),
48     combinedSecondaryVertexBJetTagsName_(Conf().getUntrackedParameter<string>
49     ("CombinedSecondaryVertexBJetTagsName","combinedSecondaryVertexBJetTags")),
50     combinedSecondaryVertexMVABJetTagsName_(Conf().getUntrackedParameter<string>
51     ("CombinedSecondaryVertexMVABJetTagsName","combinedSecondaryVertexMVABJetTags")),
52     trackCountingHighEffBJetTagsName_(Conf().getUntrackedParameter<string>
53     ("TrackCountingHighEffBJetTagsName","trackCountingHighEffBJetTags")),
54     trackCountingHighPurBJetTagsName_(Conf().getUntrackedParameter<string>
55     ("TrackCountingHighPurBJetTagsName","trackCountingHighPurBJetTags")),
56     softMuonBJetTagsName_(Conf().getUntrackedParameter<string>
57     ("SoftMuonBJetTagsName","softMuonBJetTags")),
58 bendavid 1.7 softMuonByIP3dBJetTagsName_(Conf().getUntrackedParameter<string>
59     ("SoftMuonByIP3dBJetTagsName","softMuonByIP3dBJetTags")),
60     softMuonByPtBJetTagsName_(Conf().getUntrackedParameter<string>
61     ("SoftMuonByPtBJetTagsName","softMuonByPtBJetTags")),
62     softElectronByIP3dBJetTagsName_(Conf().getUntrackedParameter<string>
63     ("SoftElectronByIP3dBJetTagsName","softElectronByIP3dBJetTags")),
64     softElectronByPtBJetTagsName_(Conf().getUntrackedParameter<string>
65     ("SoftElectronByPtBJetTagsName","softElectronByPtBJetTags")),
66 bendavid 1.1 pfCandMapName_(Conf().getUntrackedParameter<string>("pfCandMapName","pfCandMapName")),
67 bendavid 1.3 jetMapName_(Conf().getUntrackedParameter<string>("jetMapName","PFJetMap")),
68 bendavid 1.1 pfCandMap_(0),
69 bendavid 1.3 jetMap_(new mithep::PFJetMap),
70 bendavid 1.1 jets_(new mithep::PFJetArr(16))
71     {
72     // Constructor.
73     }
74    
75     //--------------------------------------------------------------------------------------------------
76     FillerPFJets::~FillerPFJets()
77     {
78     // Destructor.
79    
80     delete jets_;
81 bendavid 1.3 delete jetMap_;
82 bendavid 1.1 }
83    
84     //--------------------------------------------------------------------------------------------------
85 bendavid 1.10 void FillerPFJets::BookDataBlock(TreeWriter &tws)
86 bendavid 1.1 {
87     // Add jets branch to tree.
88    
89 loizides 1.2 tws.AddBranch(mitName_,&jets_);
90     OS()->add<mithep::PFJetArr>(jets_,mitName_);
91 bendavid 1.1
92 loizides 1.2 if (!pfCandMapName_.empty()) {
93     pfCandMap_ = OS()->get<PFCandidateMap>(pfCandMapName_);
94     if (pfCandMap_)
95     AddBranchDep(mitName_,pfCandMap_->GetBrName());
96     }
97 bendavid 1.3 if (!jetMapName_.empty()) {
98     jetMap_->SetBrName(mitName_);
99     OS()->add<PFJetMap>(jetMap_,jetMapName_);
100     }
101 bendavid 1.1 }
102    
103     //--------------------------------------------------------------------------------------------------
104     void FillerPFJets::FillDataBlock(const edm::Event &event,
105 loizides 1.2 const edm::EventSetup &setup)
106 bendavid 1.1 {
107     // Fill jets from edm collection into our collection.
108    
109     jets_->Delete();
110 bendavid 1.3 jetMap_->Reset();
111 bendavid 1.1
112 loizides 1.2 // handle for the Jet Collection
113 bendavid 1.1 Handle<reco::PFJetCollection> hJetProduct;
114     GetProduct(edmName_, hJetProduct, event);
115    
116 loizides 1.2 // handles for jet flavour matching
117 bendavid 1.1 Handle<reco::JetMatchedPartonsCollection> hPartonMatchingProduct;
118     if (flavorMatchingActive_)
119     GetProduct(flavorMatchingByReferenceName_, hPartonMatchingProduct, event);
120    
121     Handle<reco::JetTagCollection> hJetProbabilityBJetTags;
122     Handle<reco::JetTagCollection> hJetBProbabilityBJetTags;
123     Handle<reco::JetTagCollection> hSimpleSecondaryVertexBJetTags;
124     Handle<reco::JetTagCollection> hCombinedSecondaryVertexBJetTags;
125     Handle<reco::JetTagCollection> hCombinedSecondaryVertexMVABJetTags;
126     Handle<reco::JetTagCollection> hTrackCountingHighEffBJetTags;
127     Handle<reco::JetTagCollection> hTrackCountingHighPurBJetTags;
128     Handle<reco::JetTagCollection> hSoftMuonBJetTags;
129 bendavid 1.7 Handle<reco::JetTagCollection> hSoftMuonByIP3dBJetTags;
130     Handle<reco::JetTagCollection> hSoftMuonByPtBJetTags;
131     Handle<reco::JetTagCollection> hSoftElectronByIP3dBJetTags;
132     Handle<reco::JetTagCollection> hSoftElectronByPtBJetTags;
133 bendavid 1.1
134     if (bTaggingActive_) {
135     GetProduct(jetProbabilityBJetTagsName_, hJetProbabilityBJetTags, event);
136     GetProduct(jetBProbabilityBJetTagsName_, hJetBProbabilityBJetTags, event);
137     GetProduct(simpleSecondaryVertexBJetTagsName_, hSimpleSecondaryVertexBJetTags, event);
138     GetProduct(combinedSecondaryVertexBJetTagsName_, hCombinedSecondaryVertexBJetTags, event);
139 loizides 1.8 GetProduct(combinedSecondaryVertexMVABJetTagsName_, hCombinedSecondaryVertexMVABJetTags, event);
140 bendavid 1.1 GetProduct(trackCountingHighEffBJetTagsName_, hTrackCountingHighEffBJetTags, event);
141     GetProduct(trackCountingHighPurBJetTagsName_, hTrackCountingHighPurBJetTags, event);
142     GetProduct(softMuonBJetTagsName_, hSoftMuonBJetTags, event);
143 bendavid 1.7 GetProduct(softMuonByIP3dBJetTagsName_, hSoftMuonByIP3dBJetTags, event);
144     GetProduct(softMuonByPtBJetTagsName_, hSoftMuonByPtBJetTags, event);
145     GetProduct(softElectronByIP3dBJetTagsName_, hSoftElectronByIP3dBJetTags, event);
146     GetProduct(softElectronByPtBJetTagsName_, hSoftElectronByPtBJetTags, event);
147 bendavid 1.1 }
148    
149     const reco::PFJetCollection inJets = *(hJetProduct.product());
150    
151     //Handles to Jet to Vertex Association
152     Handle<std::vector<double> > JV_alpha;
153     Handle<std::vector<double> > JV_beta;
154     std::vector<double>::const_iterator it_jv_alpha;
155     std::vector<double>::const_iterator it_jv_beta;
156    
157     if (jetToVertexActive_) {
158     GetProduct(jetToVertexAlphaName_, JV_alpha, event);
159     GetProduct(jetToVertexBetaName_, JV_beta, event);
160     it_jv_alpha = JV_alpha->begin();
161     it_jv_beta = JV_beta->begin();
162     }
163    
164     //Define Jet Correction Services
165     const JetCorrector* correctorL2 = 0;
166     const JetCorrector* correctorL3 = 0;
167     if (jetCorrectionsActive_) {
168     correctorL2 = JetCorrector::getJetCorrector (L2JetCorrectorName_,setup);
169     correctorL3 = JetCorrector::getJetCorrector (L3JetCorrectorName_,setup);
170     }
171    
172     // loop through all jets
173     for (reco::PFJetCollection::const_iterator inJet = inJets.begin();
174     inJet != inJets.end(); ++inJet) {
175    
176     reco::PFJetRef jetRef(hJetProduct, inJet - inJets.begin());
177     reco::JetBaseRef jetBaseRef(jetRef);
178    
179     mithep::PFJet *jet = jets_->Allocate();
180     new (jet) mithep::PFJet(inJet->p4().x(),
181     inJet->p4().y(),
182     inJet->p4().z(),
183     inJet->p4().e());
184    
185 bendavid 1.3 //add to map
186     edm::Ptr<reco::PFJet> thePtr(hJetProduct, inJet - inJets.begin());
187     jetMap_->Add(thePtr, jet);
188 sixie 1.11
189     //fill jet moments
190     jet->SetSigmaEta(TMath::Sqrt(inJet->etaetaMoment()));
191     jet->SetSigmaPhi(TMath::Sqrt(inJet->phiphiMoment()));
192    
193    
194 bendavid 1.1 //fill pfjet-specific quantities
195     jet->SetChargedHadronEnergy(inJet->chargedHadronEnergy());
196     jet->SetNeutralHadronEnergy(inJet->neutralHadronEnergy());
197     jet->SetChargedEmEnergy(inJet->chargedEmEnergy());
198     jet->SetChargedMuEnergy(inJet->chargedMuEnergy());
199     jet->SetNeutralEmEnergy(inJet->neutralEmEnergy());
200     jet->SetChargedMultiplicity(inJet->chargedMultiplicity());
201     jet->SetNeutralMultiplicity(inJet->neutralMultiplicity());
202     jet->SetMuonMultiplicity(inJet->muonMultiplicity());
203    
204     if (jetToVertexActive_) {
205     //compute alpha and beta parameter for jets
206     jet->SetAlpha((*it_jv_alpha));
207     jet->SetBeta((*it_jv_beta));
208     }
209    
210     //Jet Corrections
211     if (jetCorrectionsActive_) {
212     double L2Scale = correctorL2->correction(inJet->p4());
213     double L3Scale = correctorL3->correction(inJet->p4()*L2Scale);
214     jet->SetL2RelativeCorrectionScale(L2Scale);
215     jet->SetL3AbsoluteCorrectionScale(L3Scale);
216     jet->EnableCorrection(mithep::PFJet::L2);
217     jet->EnableCorrection(mithep::PFJet::L3);
218     }
219    
220     if (bTaggingActive_) {
221     jet->SetJetProbabilityBJetTagsDisc((*(hJetProbabilityBJetTags.product()))[jetBaseRef]);
222     jet->SetJetBProbabilityBJetTagsDisc((*(hJetBProbabilityBJetTags.product()))[jetBaseRef]);
223 loizides 1.2 jet->SetSimpleSecondaryVertexBJetTagsDisc(
224     (*(hSimpleSecondaryVertexBJetTags.product()))[jetBaseRef]);
225     jet->SetCombinedSecondaryVertexBJetTagsDisc(
226     (*(hCombinedSecondaryVertexBJetTags.product()))[jetBaseRef]);
227     jet->SetCombinedSecondaryVertexMVABJetTagsDisc(
228     (*(hCombinedSecondaryVertexMVABJetTags.product()))[jetBaseRef]);
229     jet->SetTrackCountingHighEffBJetTagsDisc(
230     (*(hTrackCountingHighEffBJetTags.product()))[jetBaseRef]);
231     jet->SetTrackCountingHighPurBJetTagsDisc(
232     (*(hTrackCountingHighPurBJetTags.product()))[jetBaseRef]);
233 bendavid 1.1 jet->SetSoftMuonBJetTagsDisc((*(hSoftMuonBJetTags.product()))[jetBaseRef]);
234 bendavid 1.7 jet->SetSoftMuonByIP3dBJetTagsDisc((*(hSoftMuonByIP3dBJetTags.product()))[jetBaseRef]);
235     jet->SetSoftMuonByPtBJetTagsDisc((*(hSoftMuonByPtBJetTags.product()))[jetBaseRef]);
236 loizides 1.8 jet->SetSoftElectronByIP3dBJetTagsDisc((*(hSoftElectronByIP3dBJetTags.product()))[jetBaseRef]);
237 bendavid 1.7 jet->SetSoftElectronByPtBJetTagsDisc((*(hSoftElectronByPtBJetTags.product()))[jetBaseRef]);
238 bendavid 1.1 }
239    
240 loizides 1.4 // get the Monte Carlo flavour matching information
241 bendavid 1.1 if (flavorMatchingActive_) {
242     unsigned int iJet = inJet - inJets.begin();
243     const reco::JetMatchedPartonsCollection *matchedPartons = hPartonMatchingProduct.product();
244     reco::MatchedPartons matchedParton = (*matchedPartons)[matchedPartons->key(iJet)];
245     int flavorPhysDef = (matchedParton.physicsDefinitionParton().isNonnull())?
246     matchedParton.physicsDefinitionParton()->pdgId():0;
247     int flavorAlgDef = (matchedParton.algoDefinitionParton().isNonnull())?
248     matchedParton.algoDefinitionParton()->pdgId():0;
249    
250     if (flavorMatchingDefinition_ == "Algorithmic") {
251     jet->SetMatchedMCFlavor(flavorAlgDef);
252     } else if(flavorMatchingDefinition_ == "Physics") {
253     jet->SetMatchedMCFlavor(flavorPhysDef);
254     } else {
255     jet->SetMatchedMCFlavor(0);
256     }
257     }
258    
259 loizides 1.4 // add PFCandidate refs
260 bendavid 1.1 if (pfCandMap_) {
261     for (uint i=0; i<inJet->numberOfDaughters(); ++i) {
262     const reco::CandidatePtr candPtr = inJet->daughterPtr(i);
263     const reco::PFCandidatePtr pfPtr(candPtr);
264     const mithep::PFCandidate *constituent = pfCandMap_->GetMit(pfPtr);
265     jet->AddPFCand(constituent);
266     }
267     }
268     if (jetToVertexActive_) {
269     it_jv_alpha++;
270     it_jv_beta++;
271     }
272     }
273     jets_->Trim();
274     }