CMS 3D CMS Logo

GlobalMuonMatchAnalyzer.cc
Go to the documentation of this file.
1 
12 
13 // user include files
19 
22 
24 
28 
29 
32 
33 
34 #include <TH2.h>
35 
36 
38 
39 {
40  iConfig = ps;
41  //now do what ever initialization is needed
44 
45  tkAssociatorToken_ = consumes<reco::TrackToTrackingParticleAssociator>(tkAssociatorName_);
46  muAssociatorToken_ = consumes<reco::TrackToTrackingParticleAssociator>(muAssociatorName_);
47 
48 subsystemname_ = iConfig.getUntrackedParameter<std::string>("subSystemFolder", "YourSubsystem") ;
53 
56 
57  tpToken_ = consumes<edm::View<reco::Track> >(tpName_);
58  tkToken_ = consumes<edm::View<reco::Track> >(tkName_);
59  staToken_ = consumes<edm::View<reco::Track> >(staName_);
60  glbToken_ = consumes<edm::View<reco::Track> >(glbName_);
61 }
62 
63 
65 {
66 
67  // do anything here that needs to be done at desctruction time
68  // (e.g. close files, deallocate resources etc.)
69 
70 }
71 
72 
73 //
74 // member functions
75 //
76 
77 // ------------ method called to for each event ------------
78 void
80 {
81  using namespace edm;
82  using namespace reco;
83 
85  iEvent.getByToken(tpToken_,tpHandle);
86  const TrackingParticleCollection tpColl = *(tpHandle.product());
87 
89  iEvent.getByToken(glbToken_,muHandle);
90  const reco::MuonTrackLinksCollection muColl = *(muHandle.product());
91 
92  Handle<View<Track> > staHandle;
93  iEvent.getByToken(staToken_,staHandle);
94 
95  Handle<View<Track> > glbHandle;
96  iEvent.getByToken(glbToken_,glbHandle);
97 
98  Handle<View<Track> > tkHandle;
99  iEvent.getByToken(tkToken_,tkHandle);
100 
102  iEvent.getByToken(tkAssociatorToken_,tkAssociator);
103 
104  // Mu Associator
106  iEvent.getByToken(muAssociatorToken_,muAssociator);
107 
108 
109  reco::RecoToSimCollection tkrecoToSimCollection = tkAssociator->associateRecoToSim(tkHandle,tpHandle);
110  reco::SimToRecoCollection tksimToRecoCollection = tkAssociator->associateSimToReco(tkHandle,tpHandle);
111 
112  reco::RecoToSimCollection starecoToSimCollection = muAssociator->associateRecoToSim(staHandle,tpHandle);
113  reco::SimToRecoCollection stasimToRecoCollection = muAssociator->associateSimToReco(staHandle,tpHandle);
114 
115  reco::RecoToSimCollection glbrecoToSimCollection = muAssociator->associateRecoToSim(glbHandle,tpHandle);
116  reco::SimToRecoCollection glbsimToRecoCollection = muAssociator->associateSimToReco(glbHandle,tpHandle);
117 
118 
119  for (TrackingParticleCollection::size_type i=0; i<tpColl.size(); ++i){
120  TrackingParticleRef tp(tpHandle,i);
121 
122  std::vector<std::pair<RefToBase<Track>, double> > rvGlb;
123  RefToBase<Track> rGlb;
124  if(glbsimToRecoCollection.find(tp) != glbsimToRecoCollection.end()){
125  rvGlb = glbsimToRecoCollection[tp];
126  if(!rvGlb.empty()) {
127  rGlb = rvGlb.begin()->first;
128  }
129  }
130 
131  std::vector<std::pair<RefToBase<Track>, double> > rvSta;
132  RefToBase<Track> rSta;
133  if(stasimToRecoCollection.find(tp) != stasimToRecoCollection.end()){
134  rvSta = stasimToRecoCollection[tp];
135  if(!rvSta.empty()) {
136  rSta = rvSta.begin()->first;
137  }
138  }
139 
140  std::vector<std::pair<RefToBase<Track>, double> > rvTk;
141  RefToBase<Track> rTk;
142  if(tksimToRecoCollection.find(tp) != tksimToRecoCollection.end()){
143  rvTk = tksimToRecoCollection[tp];
144  if(!rvTk.empty()) {
145  rTk = rvTk.begin()->first;
146  }
147  }
148 
149  if( !rvSta.empty() && !rvTk.empty() ){
150  //should have matched
151  h_shouldMatch->Fill(rTk->eta(),rTk->pt());
152  }
153 
154  for ( reco::MuonTrackLinksCollection::const_iterator links = muHandle->begin(); links != muHandle->end(); ++links ) {
155  if( rGlb == RefToBase<Track>(links->globalTrack() ) ) {
156  if( RefToBase<Track>(links->trackerTrack() ) == rTk &&
157  RefToBase<Track>(links->standAloneTrack() ) == rSta ) {
158  //goodMatch
159  h_goodMatchSim->Fill(rGlb->eta(),rGlb->pt());
160  }
161  if ( RefToBase<Track>(links->trackerTrack() ) == rTk &&
162  RefToBase<Track>(links->standAloneTrack() ) != rSta ) {
163  //tkOnlyMatch
164  h_tkOnlySim->Fill(rGlb->eta(),rGlb->pt());
165  }
166  if ( RefToBase<Track>(links->standAloneTrack() ) == rSta &&
167  RefToBase<Track>(links->trackerTrack() ) != rTk ) {
168  //staOnlyMatch
169  h_staOnlySim->Fill(rGlb->eta(),rGlb->pt());
170  }
171  }
172  }
173 
174  }
175 
177 
178  for ( reco::MuonTrackLinksCollection::const_iterator links = muHandle->begin(); links != muHandle->end(); ++links ) {
179  RefToBase<Track> glbRef = RefToBase<Track>(links->globalTrack() );
180  RefToBase<Track> staRef = RefToBase<Track>(links->standAloneTrack() );
181  RefToBase<Track> tkRef = RefToBase<Track>(links->trackerTrack() );
182 
183  std::vector<std::pair<TrackingParticleRef, double> > tp1;
184  TrackingParticleRef tp1r;
185  if(glbrecoToSimCollection.find(glbRef) != glbrecoToSimCollection.end()){
186  tp1 = glbrecoToSimCollection[glbRef];
187  if(!tp1.empty()) {
188  tp1r = tp1.begin()->first;
189  }
190  }
191 
192  std::vector<std::pair<TrackingParticleRef, double> > tp2;
193  TrackingParticleRef tp2r;
194  if(starecoToSimCollection.find(staRef) != starecoToSimCollection.end()){
195  tp2 = starecoToSimCollection[staRef];
196  if(!tp2.empty()) {
197  tp2r = tp2.begin()->first;
198  }
199  }
200 
201  std::vector<std::pair<TrackingParticleRef, double> > tp3;
202  TrackingParticleRef tp3r;
203  if(tkrecoToSimCollection.find(tkRef) != tkrecoToSimCollection.end()){
204  tp3 = tkrecoToSimCollection[tkRef];
205  if(!tp3.empty()) {
206  tp3r = tp3.begin()->first;
207  }
208  }
209 
210 
211  if(!tp1.empty()) {
212  //was reconstructed
213  h_totReco->Fill(glbRef->eta(),glbRef->pt());
214  if(tp2r == tp3r) { // && tp1r == tp3r) {
215  //came from same TP
216  h_goodMatch->Fill(glbRef->eta(),glbRef->pt());
217  } else {
218  //mis-match
219  h_fakeMatch->Fill(glbRef->eta(),glbRef->pt());
220  }
221  }
222 
223  }
224 
225 
226 }
227 
228 
229 // ------------ method called once each job just before starting event loop ------------
230 void
232 {
233 }
235 // ------------ method called once each job just after ending the event loop ------------
239 
242 
243  if( !out.empty() && dbe_ ) dbe_->save(out);
244 }
245 
246 //void GlobalMuonMatchAnalyzer::beginRun(const edm::Run&, const edm::EventSetup& setup)
247 
249  edm::Run const & iRun,
250  edm::EventSetup const & iSetup )
251 {
252  // Tk Associator
253 
254  ibooker.cd();
255  std::string dirName="Matcher/";
256  // ibooker.setCurrentFolder("RecoMuonV/Matcher");
257  ibooker.setCurrentFolder(dirName) ;
258 
259  h_shouldMatch = ibooker.book2D("h_shouldMatch","SIM associated to Tk and Sta",50,-2.5,2.5,100,0.,500.);
260  h_goodMatchSim = ibooker.book2D("h_goodMatchSim","SIM associated to Glb Sta Tk",50,-2.5,2.5,100,0.,500.);
261  h_tkOnlySim = ibooker.book2D("h_tkOnlySim","SIM associated to Glb Tk",50,-2.5,2.5,100,0.,500.);
262  h_staOnlySim = ibooker.book2D("h_staOnlySim","SIM associated to Glb Sta",50,-2.5,2.5,100,0.,500.);
263 
264  h_totReco = ibooker.book2D("h_totReco","Total Glb Reconstructed",50,-2.5,2.5,100,0.,500.);
265  h_goodMatch = ibooker.book2D("h_goodMatch","Sta and Tk from same SIM",50,-2.5,2.5,100, 0., 500.);
266  h_fakeMatch = ibooker.book2D("h_fakeMatch","Sta and Tk not from same SIM",50,-2.5,2.5,100,0.,500.);
267 
268  h_effic = ibooker.book1D("h_effic","Efficiency vs #eta",50,-2.5,2.5);
269  h_efficPt = ibooker.book1D("h_efficPt","Efficiency vs p_{T}",100,0.,100.);
270 
271  h_fake = ibooker.book1D("h_fake","Fake fraction vs #eta",50,-2.5,2.5);
272  h_fakePt = ibooker.book1D("h_fakePt","Fake fraction vs p_{T}",100,0.,100.);
273 
274 }
275 
276 
277 
279  TH2F * h1 = recoTH2->getTH2F();
280  TH1D* reco = h1->ProjectionX();
281 
282  TH2F * h2 =simTH2->getTH2F();
283  TH1D* sim = h2 ->ProjectionX();
284 
285 
286  TH1F *hEff = (TH1F*) reco->Clone();
287 
288  hEff->Divide(sim);
289 
290  hEff->SetName("tmp_"+TString(reco->GetName()));
291 
292  // Set the error accordingly to binomial statistics
293  int nBinsEta = hEff->GetNbinsX();
294  for(int bin = 1; bin <= nBinsEta; bin++) {
295  float nSimHit = sim->GetBinContent(bin);
296  float eff = hEff->GetBinContent(bin);
297  float error = 0;
298  if(nSimHit != 0 && eff <= 1) {
299  error = sqrt(eff*(1-eff)/nSimHit);
300  }
301  hEff->SetBinError(bin, error);
302  effHist->setBinContent(bin,eff);
303  effHist->setBinError(bin,error);
304  }
305 
306 }
307 
309  TH2F * h1 = recoTH2->getTH2F();
310  TH1D* reco = h1->ProjectionY();
311 
312  TH2F * h2 = simTH2->getTH2F();
313  TH1D* sim = h2 ->ProjectionY();
314 
315 
316  TH1F *hEff = (TH1F*) reco->Clone();
317 
318  hEff->Divide(sim);
319 
320  hEff->SetName("tmp_"+TString(reco->GetName()));
321 
322  // Set the error accordingly to binomial statistics
323  int nBinsPt = hEff->GetNbinsX();
324  for(int bin = 1; bin <= nBinsPt; bin++) {
325  float nSimHit = sim->GetBinContent(bin);
326  float eff = hEff->GetBinContent(bin);
327  float error = 0;
328  if(nSimHit != 0 && eff <= 1) {
329  error = sqrt(eff*(1-eff)/nSimHit);
330  }
331  hEff->SetBinError(bin, error);
332  effHist->setBinContent(bin,eff);
333  effHist->setBinError(bin,error);
334  }
335 
336 }
337 
T getUntrackedParameter(std::string const &, T const &) const
void setBinContent(int binx, double content)
set content of bin (1-D)
edm::EDGetTokenT< reco::TrackToTrackingParticleAssociator > muAssociatorToken_
std::vector< TrackingParticle > TrackingParticleCollection
const_iterator end() const
last iterator over the map (read only)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:519
const_iterator find(const key_type &k) const
find element with specified reference key
uint16_t size_type
std::vector< MuonTrackLinks > MuonTrackLinksCollection
collection of MuonTrackLinks
Definition: MuonFwd.h:22
void Fill(long long x)
int iEvent
Definition: GenABIO.cc:230
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:651
edm::EDGetTokenT< edm::View< reco::Track > > glbToken_
T sqrt(T t)
Definition: SSEVec.h:18
void computeEfficiencyPt(MonitorElement *, MonitorElement *recoTH2, MonitorElement *simTH2)
double pt() const
track transverse momentum
Definition: TrackBase.h:621
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:118
void setBinError(int binx, double error)
set uncertainty on content of bin (1-D)
edm::EDGetTokenT< edm::View< reco::Track > > staToken_
TH2F * getTH2F() const
Definition: RunManager.h:28
bin
set the eta bin as selection string.
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:279
T const * product() const
Definition: Handle.h:81
void analyze(const edm::Event &, const edm::EventSetup &) override
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:136
reco::RecoToSimCollection associateRecoToSim(const edm::Handle< edm::View< reco::Track > > &tCH, const edm::Handle< TrackingParticleCollection > &tPCH) const
compare reco to sim the handle of reco::Track and TrackingParticle collections
GlobalMuonMatchAnalyzer(const edm::ParameterSet &)
fixed size matrix
HLT enums.
edm::EDGetTokenT< edm::View< reco::Track > > tpToken_
void endRun(const edm::Run &, const edm::EventSetup &) override
reco::SimToRecoCollection associateSimToReco(const edm::Handle< edm::View< reco::Track > > &tCH, const edm::Handle< TrackingParticleCollection > &tPCH) const
compare reco to sim the handle of reco::Track and TrackingParticle collections
void computeEfficiencyEta(MonitorElement *, MonitorElement *recoTH2, MonitorElement *simTH2)
edm::EDGetTokenT< reco::TrackToTrackingParticleAssociator > tkAssociatorToken_
Definition: Run.h:43
edm::EDGetTokenT< edm::View< reco::Track > > tkToken_