ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitProd/TreeFiller/src/FillerMuons.cc
Revision: 1.36
Committed: Mon Nov 22 16:55:51 2010 UTC (14 years, 5 months ago) by bendavid
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_017, Mit_017pre3, Mit_017pre2, Mit_017pre1
Changes since 1.35: +4 -9 lines
Log Message:
Adapt to new expected hits inner

File Contents

# User Rev Content
1 bendavid 1.36 // $Id: FillerMuons.cc,v 1.35 2010/10/30 16:22:38 pharris Exp $
2 loizides 1.1
3     #include "MitProd/TreeFiller/interface/FillerMuons.h"
4     #include "DataFormats/MuonReco/interface/Muon.h"
5 sixie 1.26 #include "DataFormats/MuonReco/interface/MuonQuality.h"
6 loizides 1.1 #include "DataFormats/MuonReco/interface/MuonFwd.h"
7 sixie 1.26 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
8 bendavid 1.15 #include "DataFormats/Common/interface/RefToPtr.h"
9 bendavid 1.27 #include "DataFormats/VertexReco/interface/VertexFwd.h"
10     #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
11     #include "TrackingTools/TransientTrack/plugins/TransientTrackBuilderESProducer.h"
12     #include "TrackingTools/IPTools/interface/IPTools.h"
13 bendavid 1.28 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexTrackCompatibilityEstimator.h"
14 loizides 1.22 #include "MitAna/DataTree/interface/Names.h"
15     #include "MitAna/DataTree/interface/MuonCol.h"
16 loizides 1.1 #include "MitAna/DataTree/interface/Track.h"
17 loizides 1.22 #include "MitProd/ObjectService/interface/ObjectService.h"
18 loizides 1.1
19     using namespace std;
20     using namespace edm;
21     using namespace mithep;
22    
23 loizides 1.4 //--------------------------------------------------------------------------------------------------
24 loizides 1.18 FillerMuons::FillerMuons(const edm::ParameterSet &cfg, const char *name, bool active) :
25     BaseFiller(cfg,name,active),
26 loizides 1.1 edmName_(Conf().getUntrackedParameter<string>("edmName","muons")),
27 bendavid 1.32 expectedHitsName_(Conf().getUntrackedParameter<string>("expectedHitsName","")),
28 loizides 1.1 mitName_(Conf().getUntrackedParameter<string>("mitName",Names::gkMuonBrn)),
29 loizides 1.9 globalTrackMapName_(Conf().getUntrackedParameter<string>("globalTrackMapName","")),
30     staTrackMapName_(Conf().getUntrackedParameter<string>("staTrackMapName","")),
31     staVtxTrackMapName_(Conf().getUntrackedParameter<string>("staVtxTrackMapName","")),
32     trackerTrackMapName_(Conf().getUntrackedParameter<string>("trackerTrackMapName","")),
33 bendavid 1.19 muonMapName_(Conf().getUntrackedParameter<string>("muonMapName","")),
34 bendavid 1.27 pvEdmName_(Conf().getUntrackedParameter<string>("pvEdmName","offlinePrimaryVertices")),
35     pvBSEdmName_(Conf().getUntrackedParameter<string>("pvEdmName","offlinePrimaryVerticesWithBS")),
36 bendavid 1.19 globalTrackMap_(0),
37     standaloneTrackMap_(0),
38     standaloneVtxTrackMap_(0),
39 loizides 1.9 trackerTrackMap_(0),
40 bendavid 1.19 muonMap_(new mithep::MuonMap),
41 pharris 1.12 muons_(new mithep::MuonArr(16))
42 loizides 1.1 {
43     // Constructor.
44     }
45    
46 loizides 1.4 //--------------------------------------------------------------------------------------------------
47 loizides 1.20 FillerMuons::~FillerMuons()
48     {
49     // Destructor.
50    
51 loizides 1.5 delete muons_;
52 bendavid 1.19 delete muonMap_;
53 loizides 1.1 }
54    
55 loizides 1.4 //--------------------------------------------------------------------------------------------------
56 bendavid 1.24 void FillerMuons::BookDataBlock(TreeWriter &tws)
57 loizides 1.20 {
58 loizides 1.9 // Add muons branch to tree and get pointers to maps.
59 loizides 1.1
60 loizides 1.20 tws.AddBranch(mitName_,&muons_);
61     OS()->add<mithep::MuonArr>(muons_,mitName_);
62 bendavid 1.19
63 loizides 1.20 if (!globalTrackMapName_.empty()) {
64     globalTrackMap_ = OS()->get<TrackMap>(globalTrackMapName_);
65     if (globalTrackMap_)
66     AddBranchDep(mitName_,globalTrackMap_->GetBrName());
67     }
68     if (!staTrackMapName_.empty()) {
69     standaloneTrackMap_ = OS()->get<TrackMap>(staTrackMapName_);
70     if (standaloneTrackMap_)
71     AddBranchDep(mitName_,standaloneTrackMap_->GetBrName());
72     }
73     if (!staVtxTrackMapName_.empty()) {
74     standaloneVtxTrackMap_ = OS()->get<TrackMap>(staVtxTrackMapName_);
75     if (standaloneVtxTrackMap_)
76     AddBranchDep(mitName_,standaloneVtxTrackMap_->GetBrName());
77     }
78     if (!trackerTrackMapName_.empty()) {
79     trackerTrackMap_ = OS()->get<TrackMap>(trackerTrackMapName_);
80     if (trackerTrackMap_)
81     AddBranchDep(mitName_,trackerTrackMap_->GetBrName());
82     }
83     if (!muonMapName_.empty()) {
84     muonMap_->SetBrName(mitName_);
85     OS()->add<MuonMap>(muonMap_,muonMapName_);
86     }
87 loizides 1.1 }
88    
89 loizides 1.4 //--------------------------------------------------------------------------------------------------
90 loizides 1.1 void FillerMuons::FillDataBlock(const edm::Event &event,
91 loizides 1.4 const edm::EventSetup &setup)
92 loizides 1.1 {
93     // Fill muon information.
94    
95 bendavid 1.16 muons_->Delete();
96 bendavid 1.19 muonMap_->Reset();
97 loizides 1.1
98 loizides 1.6 Handle<reco::MuonCollection> hMuonProduct;
99 sixie 1.10 GetProduct(edmName_, hMuonProduct, event);
100 loizides 1.6 const reco::MuonCollection inMuons = *(hMuonProduct.product());
101 loizides 1.1
102 bendavid 1.27 edm::Handle<reco::VertexCollection> hVertex;
103     event.getByLabel(pvEdmName_, hVertex);
104     const reco::VertexCollection *pvCol = hVertex.product();
105    
106     edm::Handle<reco::VertexCollection> hVertexBS;
107     event.getByLabel(pvBSEdmName_, hVertexBS);
108     const reco::VertexCollection *pvBSCol = hVertexBS.product();
109 bendavid 1.32
110 bendavid 1.27 edm::ESHandle<TransientTrackBuilder> hTransientTrackBuilder;
111     setup.get<TransientTrackRecord>().get("TransientTrackBuilder",hTransientTrackBuilder);
112     const TransientTrackBuilder *transientTrackBuilder = hTransientTrackBuilder.product();
113    
114 bendavid 1.28 KalmanVertexTrackCompatibilityEstimator<5> kalmanEstimator;
115    
116 sixie 1.10 for (reco::MuonCollection::const_iterator iM = inMuons.begin(); iM != inMuons.end(); ++iM) {
117 bendavid 1.25 mithep::Muon* outMuon = muons_->AddNew();
118    
119     outMuon->SetPtEtaPhi(iM->pt(),iM->eta(),iM->phi());
120     outMuon->SetCharge(iM->charge());
121    
122 sixie 1.10 outMuon->SetIsoR03SumPt(iM->isolationR03().sumPt);
123     outMuon->SetIsoR03EmEt(iM->isolationR03().emEt);
124     outMuon->SetIsoR03HadEt(iM->isolationR03().hadEt);
125     outMuon->SetIsoR03HoEt(iM->isolationR03().hoEt);
126     outMuon->SetIsoR03NTracks(iM->isolationR03().nTracks);
127     outMuon->SetIsoR03NJets(iM->isolationR03().nJets);
128     outMuon->SetIsoR05SumPt(iM->isolationR05().sumPt);
129     outMuon->SetIsoR05EmEt(iM->isolationR05().emEt);
130     outMuon->SetIsoR05HadEt(iM->isolationR05().hadEt);
131     outMuon->SetIsoR05HoEt(iM->isolationR05().hoEt);
132     outMuon->SetIsoR05NTracks(iM->isolationR05().nTracks);
133     outMuon->SetIsoR05NJets(iM->isolationR05().nJets);
134     outMuon->SetEmEnergy(iM->calEnergy().em);
135     outMuon->SetHadEnergy(iM->calEnergy().had);
136     outMuon->SetHoEnergy(iM->calEnergy().ho);
137     outMuon->SetEmS9Energy(iM->calEnergy().emS9);
138     outMuon->SetHadS9Energy(iM->calEnergy().hadS9);
139     outMuon->SetHoS9Energy(iM->calEnergy().hoS9);
140 bendavid 1.17 outMuon->SetIsGlobalMuon(iM->isGlobalMuon());
141     outMuon->SetIsTrackerMuon(iM->isTrackerMuon());
142     outMuon->SetIsStandaloneMuon(iM->isStandAloneMuon());
143     outMuon->SetIsCaloMuon(iM->isCaloMuon());
144 pharris 1.11
145 sixie 1.26 //save results from the muon selector in the MuonQuality bitmask
146     outMuon->Quality().SetQuality(MuonQuality::All);
147     if (muon::isGoodMuon(*iM,muon::AllGlobalMuons))
148     outMuon->Quality().SetQuality(MuonQuality::AllGlobalMuons);
149     if (muon::isGoodMuon(*iM,muon::AllStandAloneMuons))
150     outMuon->Quality().SetQuality(MuonQuality::AllStandAloneMuons);
151     if (muon::isGoodMuon(*iM,muon::AllTrackerMuons))
152     outMuon->Quality().SetQuality(MuonQuality::AllTrackerMuons);
153     if (muon::isGoodMuon(*iM,muon::TrackerMuonArbitrated))
154     outMuon->Quality().SetQuality(MuonQuality::TrackerMuonArbitrated);
155     if (muon::isGoodMuon(*iM,muon::AllArbitrated))
156     outMuon->Quality().SetQuality(MuonQuality::AllArbitrated);
157     if (muon::isGoodMuon(*iM,muon::GlobalMuonPromptTight))
158     outMuon->Quality().SetQuality(MuonQuality::GlobalMuonPromptTight);
159     if (muon::isGoodMuon(*iM,muon::TMLastStationLoose))
160     outMuon->Quality().SetQuality(MuonQuality::TMLastStationLoose);
161     if (muon::isGoodMuon(*iM,muon::TMLastStationTight))
162     outMuon->Quality().SetQuality(MuonQuality::TMLastStationTight);
163     if (muon::isGoodMuon(*iM,muon::TM2DCompatibilityLoose))
164     outMuon->Quality().SetQuality(MuonQuality::TM2DCompatibilityLoose);
165     if (muon::isGoodMuon(*iM,muon::TM2DCompatibilityTight))
166     outMuon->Quality().SetQuality(MuonQuality::TM2DCompatibilityTight);
167     if (muon::isGoodMuon(*iM,muon::TMOneStationLoose))
168     outMuon->Quality().SetQuality(MuonQuality::TMOneStationLoose);
169     if (muon::isGoodMuon(*iM,muon::TMOneStationTight))
170     outMuon->Quality().SetQuality(MuonQuality::TMOneStationTight);
171     if (muon::isGoodMuon(*iM,muon::TMLastStationOptimizedLowPtLoose))
172     outMuon->Quality().SetQuality(MuonQuality::TMLastStationOptimizedLowPtLoose);
173     if (muon::isGoodMuon(*iM,muon::TMLastStationOptimizedLowPtTight))
174     outMuon->Quality().SetQuality(MuonQuality::TMLastStationOptimizedLowPtTight);
175     if (muon::isGoodMuon(*iM,muon::GMTkChiCompatibility))
176     outMuon->Quality().SetQuality(MuonQuality::GMTkChiCompatibility);
177     if (muon::isGoodMuon(*iM,muon::GMStaChiCompatibility))
178     outMuon->Quality().SetQuality(MuonQuality::GMStaChiCompatibility);
179     if (muon::isGoodMuon(*iM,muon::GMTkKinkTight))
180     outMuon->Quality().SetQuality(MuonQuality::GMTkKinkTight);
181     if (muon::isGoodMuon(*iM,muon::TMLastStationAngLoose))
182     outMuon->Quality().SetQuality(MuonQuality::TMLastStationAngLoose);
183     if (muon::isGoodMuon(*iM,muon::TMLastStationAngTight))
184     outMuon->Quality().SetQuality(MuonQuality::TMLastStationAngTight);
185     if (muon::isGoodMuon(*iM,muon::TMOneStationAngLoose))
186     outMuon->Quality().SetQuality(MuonQuality::TMOneStationAngLoose);
187     if (muon::isGoodMuon(*iM,muon::TMOneStationAngTight))
188     outMuon->Quality().SetQuality(MuonQuality::TMOneStationAngTight);
189     if (muon::isGoodMuon(*iM,muon::TMLastStationOptimizedBarrelLowPtLoose))
190     outMuon->Quality().SetQuality(MuonQuality::TMLastStationOptimizedBarrelLowPtLoose);
191     if (muon::isGoodMuon(*iM,muon::TMLastStationOptimizedBarrelLowPtTight))
192     outMuon->Quality().SetQuality(MuonQuality::TMLastStationOptimizedBarrelLowPtTight);
193    
194    
195 bendavid 1.30 if (globalTrackMap_ && iM->combinedMuon().isNonnull()) {
196 bendavid 1.15 outMuon->SetGlobalTrk(globalTrackMap_->GetMit(refToPtr(iM->combinedMuon())));
197 bendavid 1.30 outMuon->SetNValidHits (iM->globalTrack()->hitPattern().numberOfValidMuonHits());
198     }
199 loizides 1.1 if (standaloneTrackMap_ && standaloneVtxTrackMap_ && iM->standAloneMuon().isNonnull()) {
200     Int_t refProductId = iM->standAloneMuon().id().id();
201     if ( refProductId == standaloneVtxTrackMap_->GetEdmProductId())
202 bendavid 1.15 outMuon->SetStandaloneTrk(standaloneVtxTrackMap_->GetMit(refToPtr(iM->standAloneMuon())));
203 loizides 1.1 else if ( refProductId == standaloneTrackMap_->GetEdmProductId())
204 bendavid 1.15 outMuon->SetStandaloneTrk(standaloneTrackMap_->GetMit(refToPtr(iM->standAloneMuon())));
205 loizides 1.1 else throw edm::Exception(edm::errors::Configuration, "FillerMuons:FillDataBlock()\n")
206 loizides 1.14 << "Error! Track reference in unmapped collection " << edmName_ << endl;
207 loizides 1.1 }
208 bendavid 1.30 if (trackerTrackMap_ && iM->track().isNonnull()) {
209 bendavid 1.15 outMuon->SetTrackerTrk(trackerTrackMap_->GetMit(refToPtr(iM->track())));
210 bendavid 1.30 }
211 loizides 1.14
212 bendavid 1.27 //compute impact parameter with respect to PV
213     if (iM->track().isNonnull()) {
214     const reco::TransientTrack &tt = transientTrackBuilder->build(iM->track());
215    
216     const std::pair<bool,Measurement1D> &d0pv = IPTools::absoluteTransverseImpactParameter(tt,pvCol->at(0));
217     if (d0pv.first) {
218     outMuon->SetD0PV(d0pv.second.value());
219     outMuon->SetD0PVErr(d0pv.second.error());
220     }
221 bendavid 1.28 else {
222     outMuon->SetD0PV(-99.0);
223     }
224 bendavid 1.27
225     const std::pair<bool,Measurement1D> &ip3dpv = IPTools::absoluteImpactParameter3D(tt,pvCol->at(0));
226     if (ip3dpv.first) {
227     outMuon->SetIp3dPV(ip3dpv.second.value());
228     outMuon->SetIp3dPVErr(ip3dpv.second.error());
229     }
230 bendavid 1.28 else {
231     outMuon->SetIp3dPV(-99.0);
232     }
233 bendavid 1.27
234     const std::pair<bool,Measurement1D> &d0pvbs = IPTools::absoluteTransverseImpactParameter(tt,pvBSCol->at(0));
235     if (d0pvbs.first) {
236     outMuon->SetD0PVBS(d0pvbs.second.value());
237     outMuon->SetD0PVBSErr(d0pvbs.second.error());
238     }
239 bendavid 1.28 else {
240     outMuon->SetD0PVBS(-99.0);
241     }
242 bendavid 1.27
243     const std::pair<bool,Measurement1D> &ip3dpvbs = IPTools::absoluteImpactParameter3D(tt,pvBSCol->at(0));
244     if (ip3dpvbs.first) {
245     outMuon->SetIp3dPVBS(ip3dpvbs.second.value());
246     outMuon->SetIp3dPVBSErr(ip3dpvbs.second.error());
247     }
248 bendavid 1.28 else {
249     outMuon->SetIp3dPVBS(-99.0);
250     }
251    
252     //compute compatibility with PV using taking into account also the case where muon track
253     //was included in the vertex fit
254     if (iM->track()->extra().isAvailable()) {
255    
256     const std::pair<bool,double> &pvCompat = kalmanEstimator.estimate(pvCol->at(0),tt);
257     if (pvCompat.first) {
258     outMuon->SetPVCompatibility(pvCompat.second);
259     }
260     else {
261     outMuon->SetPVCompatibility(-99.0);
262     }
263    
264     const std::pair<bool,double> &pvbsCompat = kalmanEstimator.estimate(pvBSCol->at(0),tt);
265     if (pvbsCompat.first) {
266     outMuon->SetPVBSCompatibility(pvbsCompat.second);
267     }
268     else {
269     outMuon->SetPVBSCompatibility(-99.0);
270     }
271     }
272 bendavid 1.27
273     }
274 bendavid 1.31
275 loizides 1.14 outMuon->SetNChambers (iM->numberOfChambers());
276 pharris 1.11 outMuon->SetStationMask(iM->stationMask(reco::Muon::SegmentAndTrackArbitration));
277 ceballos 1.33 outMuon->SetNMatches (iM->numberOfMatches());
278 loizides 1.14 for(int i0 = 0; i0 < 4; i0++) {
279 loizides 1.20 // DTs
280 pharris 1.11 outMuon->SetDX(i0, iM->dX(i0+1,1));
281     outMuon->SetDY(i0, iM->dY(i0+1,1));
282     outMuon->SetPullX(i0, iM->pullX(i0+1,1));
283     outMuon->SetPullY(i0, iM->pullY(i0+1,1));
284     outMuon->SetTrackDist(i0, iM->trackDist(i0+1,1));
285     outMuon->SetTrackDistErr(i0, iM->trackDistErr(i0+1,1));
286 pharris 1.34 outMuon->SetNSegments(i0, NumberOfSegments(&(*iM),i0+1,1));
287 loizides 1.20 // CSCs
288 pharris 1.11 outMuon->SetDX(4+i0, iM->dX (i0+1,2));
289     outMuon->SetDY(4+i0, iM->dY (i0+1,2));
290     outMuon->SetPullX(4+i0, iM->pullX (i0+1,2));
291     outMuon->SetPullY(4+i0, iM->pullY (i0+1,2));
292     outMuon->SetTrackDist(4+i0, iM->trackDist(i0+1,2));
293     outMuon->SetTrackDistErr(4+i0,iM->trackDistErr(i0+1,2));
294 pharris 1.34 outMuon->SetNSegments(4+i0, NumberOfSegments(&(*iM),i0+1,2));
295 pharris 1.11 }
296 bendavid 1.19
297 bendavid 1.32 reco::MuonRef theRef(hMuonProduct, iM - inMuons.begin());
298     // fill corrected expected inner hits
299 bendavid 1.36 if (iM->innerTrack().isNonnull()) {
300     outMuon->SetCorrectedNExpectedHitsInner(iM->innerTrack()->trackerExpectedHitsInner().numberOfHits());
301 bendavid 1.32 }
302 bendavid 1.36
303 loizides 1.20 // add muon to map
304 bendavid 1.19 edm::Ptr<reco::Muon> thePtr(hMuonProduct, iM - inMuons.begin());
305     muonMap_->Add(thePtr, outMuon);
306 bendavid 1.21
307     if (verbose_>1) {
308     if (!outMuon->HasGlobalTrk() && !outMuon->HasStandaloneTrk()) {
309     printf("mithep::Muon, pt=%5f, eta=%5f, phi=%5f, mass=%5f\n",outMuon->Pt(),outMuon->Eta(),outMuon->Phi(), outMuon->Mass());
310     printf(" reco::Muon, pt=%5f, eta=%5f, phi=%5f, mass=%5f\n",iM->pt(),iM->eta(),iM->phi(),iM->mass());
311     }
312     }
313 bendavid 1.19
314 pharris 1.11 }
315 loizides 1.1 muons_->Trim();
316     }
317 pharris 1.34
318     int FillerMuons::NumberOfSegments(const reco::Muon *iM, int station, int muonSubdetId, reco::Muon::ArbitrationType type ) {
319     if(!iM->isMatchesValid()) return 0;
320     int segments(0);
321    
322     for( std::vector<reco::MuonChamberMatch>::const_iterator chamberMatch = iM->matches().begin();
323     chamberMatch != iM->matches().end(); chamberMatch++ )
324     {
325     if(chamberMatch->segmentMatches.empty()) continue;
326     if(!(chamberMatch->station()==station && chamberMatch->detector()==muonSubdetId)) continue;
327     if(type == reco::Muon::NoArbitration) {
328     segments += chamberMatch->segmentMatches.size();
329     continue;
330     }
331     for( std::vector<reco::MuonSegmentMatch>::const_iterator segmentMatch = chamberMatch->segmentMatches.begin();
332     segmentMatch != chamberMatch->segmentMatches.end(); segmentMatch++ )
333     {
334     if(type == reco::Muon::SegmentArbitration)
335 pharris 1.35 if(segmentMatch->isMask(reco::MuonSegmentMatch::BestInChamberByDR)) {
336 pharris 1.34 segments++;
337     break;
338     }
339     if(type == reco::Muon::SegmentAndTrackArbitration)
340 pharris 1.35 if(segmentMatch->isMask(reco::MuonSegmentMatch::BestInChamberByDR) &&
341 pharris 1.34 segmentMatch->isMask(reco::MuonSegmentMatch::BelongsToTrackByDR)) {
342     segments++;
343     break;
344     }
345     if(type == reco::Muon::SegmentAndTrackArbitrationCleaned)
346 pharris 1.35 if(segmentMatch->isMask(reco::MuonSegmentMatch::BestInChamberByDR) &&
347 pharris 1.34 segmentMatch->isMask(reco::MuonSegmentMatch::BelongsToTrackByDR) &&
348     segmentMatch->isMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning)) {
349     segments++;
350     break;
351     }
352     if(type > 1<<7)
353     if(segmentMatch->isMask(type)) {
354     segments++;
355     break;
356     }
357     }
358     }
359     return segments;
360     }
361