CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MuonTrackValidator.cc
Go to the documentation of this file.
3 
6 
19 
20 #include "TMath.h"
21 #include <TF1.h>
22 
23 using namespace std;
24 using namespace edm;
25 
26 
28 
29 
30  int j=0;
31  for (unsigned int ww=0;ww<associators.size();ww++){
32  for (unsigned int www=0;www<label.size();www++){
33 
34  ibooker.cd();
35  InputTag algo = label[www];
36  string dirName=dirName_;
37  if (algo.process()!="")
38  dirName+=algo.process()+"_";
39  if(algo.label()!="")
40  dirName+=algo.label();
41  if(algo.instance()!="")
42  dirName+=("_"+algo.instance());
43  if (dirName.find("Tracks")<dirName.length()){
44  dirName.replace(dirName.find("Tracks"),6,"Trks");
45  }
46  if (dirName.find("UpdatedAtVtx")<dirName.length()){
47  dirName.replace(dirName.find("UpdatedAtVtx"),12,"UpdAtVtx");
48  }
49  string assoc= associators[ww];
50  if (assoc.find("tpToTkmuTrackAssociation")<assoc.length()){
51  dirName+="_TkAsso";
52  }
53  std::replace(dirName.begin(), dirName.end(), ':', '_');
54  ibooker.setCurrentFolder(dirName.c_str());
55 
56  setUpVectors();
57 
58  ibooker.goUp();
59  string subDirName = dirName + "/simulation";
60  ibooker.setCurrentFolder(subDirName.c_str());
61  h_ptSIM.push_back( ibooker.book1D("ptSIM", "generated p_{t}", 5500, 0, 110 ) );
62  h_etaSIM.push_back( ibooker.book1D("etaSIM", "generated pseudorapidity", 500, -2.5, 2.5 ) );
63  h_tracksSIM.push_back( ibooker.book1D("tracksSIM","number of simulated tracks",200,-0.5,99.5) );
64  h_vertposSIM.push_back( ibooker.book1D("vertposSIM","Transverse position of sim vertices",100,0.,120.) );
65 
66  ibooker.cd();
67  ibooker.setCurrentFolder(dirName.c_str());
68  h_tracks.push_back( ibooker.book1D("tracks","number of reconstructed tracks",200,-0.5,19.5) );
69  h_fakes.push_back( ibooker.book1D("fakes","number of fake reco tracks",20,-0.5,19.5) );
70  h_charge.push_back( ibooker.book1D("charge","charge",3,-1.5,1.5) );
71  h_hits.push_back( ibooker.book1D("hits", "number of hits per track", nintHit,minHit,maxHit ) );
72  h_losthits.push_back( ibooker.book1D("losthits", "number of lost hits per track", nintHit,minHit,maxHit) );
73  h_nchi2.push_back( ibooker.book1D("chi2", "normalized #chi^{2}", 200, 0, 20 ) );
74  h_nchi2_prob.push_back( ibooker.book1D("chi2_prob", "normalized #chi^{2} probability",100,0,1));
75 
77  h_recoeta.push_back( ibooker.book1D("num_reco_eta","N of reco track vs eta",nint,min,max) );
78  h_assoceta.push_back( ibooker.book1D("num_assoc(simToReco)_eta","N of associated tracks (simToReco) vs eta",nint,min,max) );
79  h_assoc2eta.push_back( ibooker.book1D("num_assoc(recoToSim)_eta","N of associated (recoToSim) tracks vs eta",nint,min,max) );
80  h_simuleta.push_back( ibooker.book1D("num_simul_eta","N of simulated tracks vs eta",nint,min,max) );
81  h_recopT.push_back( ibooker.book1D("num_reco_pT","N of reco track vs pT",nintpT,minpT,maxpT) );
82  h_assocpT.push_back( ibooker.book1D("num_assoc(simToReco)_pT","N of associated tracks (simToReco) vs pT",nintpT,minpT,maxpT) );
83  h_assoc2pT.push_back( ibooker.book1D("num_assoc(recoToSim)_pT","N of associated (recoToSim) tracks vs pT",nintpT,minpT,maxpT) );
84  h_simulpT.push_back( ibooker.book1D("num_simul_pT","N of simulated tracks vs pT",nintpT,minpT,maxpT) );
85  //
86  h_recohit.push_back( ibooker.book1D("num_reco_hit","N of reco track vs hit",nintHit,minHit,maxHit) );
87  h_assochit.push_back( ibooker.book1D("num_assoc(simToReco)_hit","N of associated tracks (simToReco) vs hit",nintHit,minHit,maxHit) );
88  h_assoc2hit.push_back( ibooker.book1D("num_assoc(recoToSim)_hit","N of associated (recoToSim) tracks vs hit",nintHit,minHit,maxHit) );
89  h_simulhit.push_back( ibooker.book1D("num_simul_hit","N of simulated tracks vs hit",nintHit,minHit,maxHit) );
90  //
91  h_recophi.push_back( ibooker.book1D("num_reco_phi","N of reco track vs phi",nintPhi,minPhi,maxPhi) );
92  h_assocphi.push_back( ibooker.book1D("num_assoc(simToReco)_phi","N of associated tracks (simToReco) vs phi",nintPhi,minPhi,maxPhi) );
93  h_assoc2phi.push_back( ibooker.book1D("num_assoc(recoToSim)_phi","N of associated (recoToSim) tracks vs phi",nintPhi,minPhi,maxPhi) );
94  h_simulphi.push_back( ibooker.book1D("num_simul_phi","N of simulated tracks vs phi",nintPhi,minPhi,maxPhi) );
95 
96  h_recodxy.push_back( ibooker.book1D("num_reco_dxy","N of reco track vs dxy",nintDxy,minDxy,maxDxy) );
97  h_assocdxy.push_back( ibooker.book1D("num_assoc(simToReco)_dxy","N of associated tracks (simToReco) vs dxy",nintDxy,minDxy,maxDxy) );
98  h_assoc2dxy.push_back( ibooker.book1D("num_assoc(recoToSim)_dxy","N of associated (recoToSim) tracks vs dxy",nintDxy,minDxy,maxDxy) );
99  h_simuldxy.push_back( ibooker.book1D("num_simul_dxy","N of simulated tracks vs dxy",nintDxy,minDxy,maxDxy) );
100 
101  h_recodz.push_back( ibooker.book1D("num_reco_dz","N of reco track vs dz",nintDz,minDz,maxDz) );
102  h_assocdz.push_back( ibooker.book1D("num_assoc(simToReco)_dz","N of associated tracks (simToReco) vs dz",nintDz,minDz,maxDz) );
103  h_assoc2dz.push_back( ibooker.book1D("num_assoc(recoToSim)_dz","N of associated (recoToSim) tracks vs dz",nintDz,minDz,maxDz) );
104  h_simuldz.push_back( ibooker.book1D("num_simul_dz","N of simulated tracks vs dz",nintDz,minDz,maxDz) );
105 
106  h_assocvertpos.push_back( ibooker.book1D("num_assoc(simToReco)_vertpos","N of associated tracks (simToReco) vs transverse vert position",nintVertpos,minVertpos,maxVertpos) );
107  h_simulvertpos.push_back( ibooker.book1D("num_simul_vertpos","N of simulated tracks vs transverse vert position",nintVertpos,minVertpos,maxVertpos) );
108 
109  h_assoczpos.push_back( ibooker.book1D("num_assoc(simToReco)_zpos","N of associated tracks (simToReco) vs z vert position",nintZpos,minZpos,maxZpos) );
110  h_simulzpos.push_back( ibooker.book1D("num_simul_zpos","N of simulated tracks vs z vert position",nintZpos,minZpos,maxZpos) );
111 
112 
114 
115  h_eta.push_back( ibooker.book1D("eta", "pseudorapidity residue", 1000, -0.1, 0.1 ) );
116  h_pt.push_back( ibooker.book1D("pullPt", "pull of p_{t}", 100, -10, 10 ) );
117  h_pullTheta.push_back( ibooker.book1D("pullTheta","pull of #theta parameter",250,-25,25) );
118  h_pullPhi.push_back( ibooker.book1D("pullPhi","pull of #phi parameter",250,-25,25) );
119  h_pullDxy.push_back( ibooker.book1D("pullDxy","pull of dxy parameter",250,-25,25) );
120  h_pullDz.push_back( ibooker.book1D("pullDz","pull of dz parameter",250,-25,25) );
121  h_pullQoverp.push_back( ibooker.book1D("pullQoverp","pull of qoverp parameter",250,-25,25) );
122 
123  if (associators[ww]=="trackAssociatorByChi2"){
124  h_assochi2.push_back( ibooker.book1D("assocChi2","track association #chi^{2}",1000000,0,100000) );
125  h_assochi2_prob.push_back(ibooker.book1D("assocChi2_prob","probability of association #chi^{2}",100,0,1));
126  } else if (associators[ww]=="trackAssociatorByHits"){
127  h_assocFraction.push_back( ibooker.book1D("assocFraction","fraction of shared hits",200,0,2) );
128  h_assocSharedHit.push_back(ibooker.book1D("assocSharedHit","number of shared hits",20,0,20));
129  }
130 
131  chi2_vs_nhits.push_back( ibooker.book2D("chi2_vs_nhits","#chi^{2} vs nhits",25,0,25,100,0,10) );
132  h_chi2meanhitsh.push_back( ibooker.bookProfile("chi2mean_vs_nhits","mean #chi^{2} vs nhits",25,0,25,100,0,10) );
133 
134  etares_vs_eta.push_back( ibooker.book2D("etares_vs_eta","etaresidue vs eta",nint,min,max,200,-0.1,0.1) );
135  nrec_vs_nsim.push_back( ibooker.book2D("nrec_vs_nsim","nrec vs nsim",20,-0.5,19.5,20,-0.5,19.5) );
136 
137  chi2_vs_eta.push_back( ibooker.book2D("chi2_vs_eta","chi2_vs_eta",nint,min,max, 200, 0, 20 ));
138  h_chi2meanh.push_back( ibooker.bookProfile("chi2mean","mean #chi^{2} vs #eta",nint,min,max, 200, 0, 20) );
139  chi2_vs_phi.push_back( ibooker.book2D("chi2_vs_phi","#chi^{2} vs #phi",nintPhi,minPhi,maxPhi, 200, 0, 20 ) );
140  h_chi2mean_vs_phi.push_back( ibooker.bookProfile("chi2mean_vs_phi","mean of #chi^{2} vs #phi",nintPhi,minPhi,maxPhi, 200, 0, 20) );
141 
142  nhits_vs_eta.push_back( ibooker.book2D("nhits_vs_eta","nhits vs eta",nint,min,max,nintHit,minHit,maxHit) );
143  nDThits_vs_eta.push_back( ibooker.book2D("nDThits_vs_eta","# DT hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
144  nCSChits_vs_eta.push_back( ibooker.book2D("nCSChits_vs_eta","# CSC hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
145  nRPChits_vs_eta.push_back( ibooker.book2D("nRPChits_vs_eta","# RPC hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
146 
147  h_DThits_eta.push_back( ibooker.bookProfile("DThits_eta","mean # DT hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
148  h_CSChits_eta.push_back( ibooker.bookProfile("CSChits_eta","mean # CSC hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
149  h_RPChits_eta.push_back( ibooker.bookProfile("RPChits_eta","mean # RPC hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
150  h_hits_eta.push_back( ibooker.bookProfile("hits_eta","mean #hits vs eta",nint,min,max,nintHit,minHit,maxHit) );
151  nhits_vs_phi.push_back( ibooker.book2D("nhits_vs_phi","#hits vs #phi",nintPhi,minPhi,maxPhi,nintHit,minHit,maxHit) );
152  h_hits_phi.push_back( ibooker.bookProfile("hits_phi","mean #hits vs #phi",nintPhi,minPhi,maxPhi, nintHit,minHit,maxHit) );
153 
154  nlosthits_vs_eta.push_back( ibooker.book2D("nlosthits_vs_eta","nlosthits vs eta",nint,min,max,nintHit,minHit,maxHit) );
155  h_losthits_eta.push_back( ibooker.bookProfile("losthits_eta","losthits_eta",nint,min,max,nintHit,minHit,maxHit) );
156 
157  ptres_vs_eta.push_back(ibooker.book2D("ptres_vs_eta","ptres_vs_eta",nint,min,max, ptRes_nbin, ptRes_rangeMin, ptRes_rangeMax));
158  ptres_vs_phi.push_back( ibooker.book2D("ptres_vs_phi","p_{t} res vs #phi",nintPhi,minPhi,maxPhi, ptRes_nbin, ptRes_rangeMin, ptRes_rangeMax));
159  ptres_vs_pt.push_back(ibooker.book2D("ptres_vs_pt","ptres_vs_pt",nintpT,minpT,maxpT, ptRes_nbin, ptRes_rangeMin, ptRes_rangeMax));
160 
161  cotThetares_vs_eta.push_back(ibooker.book2D("cotThetares_vs_eta","cotThetares_vs_eta",nint,min,max,cotThetaRes_nbin, cotThetaRes_rangeMin, cotThetaRes_rangeMax));
162  cotThetares_vs_pt.push_back(ibooker.book2D("cotThetares_vs_pt","cotThetares_vs_pt",nintpT,minpT,maxpT, cotThetaRes_nbin, cotThetaRes_rangeMin, cotThetaRes_rangeMax));
163 
164  phires_vs_eta.push_back(ibooker.book2D("phires_vs_eta","phires_vs_eta",nint,min,max, phiRes_nbin, phiRes_rangeMin, phiRes_rangeMax));
165  phires_vs_pt.push_back(ibooker.book2D("phires_vs_pt","phires_vs_pt",nintpT,minpT,maxpT, phiRes_nbin, phiRes_rangeMin, phiRes_rangeMax));
166  phires_vs_phi.push_back(ibooker.book2D("phires_vs_phi","#phi res vs #phi",nintPhi,minPhi,maxPhi,phiRes_nbin, phiRes_rangeMin, phiRes_rangeMax));
167 
168  dxyres_vs_eta.push_back(ibooker.book2D("dxyres_vs_eta","dxyres_vs_eta",nint,min,max,dxyRes_nbin, dxyRes_rangeMin, dxyRes_rangeMax));
169  dxyres_vs_pt.push_back( ibooker.book2D("dxyres_vs_pt","dxyres_vs_pt",nintpT,minpT,maxpT,dxyRes_nbin, dxyRes_rangeMin, dxyRes_rangeMax));
170 
171  dzres_vs_eta.push_back(ibooker.book2D("dzres_vs_eta","dzres_vs_eta",nint,min,max,dzRes_nbin, dzRes_rangeMin, dzRes_rangeMax));
172  dzres_vs_pt.push_back(ibooker.book2D("dzres_vs_pt","dzres_vs_pt",nintpT,minpT,maxpT,dzRes_nbin, dzRes_rangeMin, dzRes_rangeMax));
173 
174  ptmean_vs_eta_phi.push_back(ibooker.bookProfile2D("ptmean_vs_eta_phi","mean p_{t} vs #eta and #phi",nintPhi,minPhi,maxPhi,nint,min,max,1000,0,1000));
175  phimean_vs_eta_phi.push_back(ibooker.bookProfile2D("phimean_vs_eta_phi","mean #phi vs #eta and #phi",nintPhi,minPhi,maxPhi,nint,min,max,nintPhi,minPhi,maxPhi));
176 
177  //pulls of track params vs eta: to be used with fitslicesytool
178  dxypull_vs_eta.push_back(ibooker.book2D("dxypull_vs_eta","dxypull_vs_eta",nint,min,max,100,-10,10));
179  ptpull_vs_eta.push_back(ibooker.book2D("ptpull_vs_eta","ptpull_vs_eta",nint,min,max,100,-10,10));
180  dzpull_vs_eta.push_back(ibooker.book2D("dzpull_vs_eta","dzpull_vs_eta",nint,min,max,100,-10,10));
181  phipull_vs_eta.push_back(ibooker.book2D("phipull_vs_eta","phipull_vs_eta",nint,min,max,100,-10,10));
182  thetapull_vs_eta.push_back(ibooker.book2D("thetapull_vs_eta","thetapull_vs_eta",nint,min,max,100,-10,10));
183 
184  //pulls of track params vs phi
185  ptpull_vs_phi.push_back(ibooker.book2D("ptpull_vs_phi","p_{t} pull vs #phi",nintPhi,minPhi,maxPhi,100,-10,10));
186  phipull_vs_phi.push_back(ibooker.book2D("phipull_vs_phi","#phi pull vs #phi",nintPhi,minPhi,maxPhi,100,-10,10));
187  thetapull_vs_phi.push_back(ibooker.book2D("thetapull_vs_phi","#theta pull vs #phi",nintPhi,minPhi,maxPhi,100,-10,10));
188 
189  nrecHit_vs_nsimHit_sim2rec.push_back( ibooker.book2D("nrecHit_vs_nsimHit_sim2rec","nrecHit vs nsimHit (Sim2RecAssoc)",nintHit,minHit,maxHit, nintHit,minHit,maxHit ));
190  nrecHit_vs_nsimHit_rec2sim.push_back( ibooker.book2D("nrecHit_vs_nsimHit_rec2sim","nrecHit vs nsimHit (Rec2simAssoc)",nintHit,minHit,maxHit, nintHit,minHit,maxHit ));
191 
192  if (MABH) {
193 h_PurityVsQuality.push_back( ibooker.book2D("PurityVsQuality","Purity vs Quality (MABH)",20,0.01,1.01,20,0.01,1.01) );
194 h_assoceta_Quality05.push_back( ibooker.book1D("num_assoc(simToReco)_eta_Q05","N of associated tracks (simToReco) vs eta (Quality>0.5)",nint,min,max) );
195 h_assoceta_Quality075.push_back( ibooker.book1D("num_assoc(simToReco)_eta_Q075","N of associated tracks (simToReco) vs eta (Quality>0.75)",nint,min,max) );
196 h_assocpT_Quality05.push_back( ibooker.book1D("num_assoc(simToReco)_pT_Q05","N of associated tracks (simToReco) vs pT (Quality>0.5)",nintpT,minpT,maxpT) );
197 h_assocpT_Quality075.push_back( ibooker.book1D("num_assoc(simToReco)_pT_Q075","N of associated tracks (simToReco) vs pT (Quality>0.75)",nintpT,minpT,maxpT) );
198 h_assocphi_Quality05.push_back( ibooker.book1D("num_assoc(simToReco)_phi_Q05","N of associated tracks (simToReco) vs phi (Quality>0.5)",nintPhi,minPhi,maxPhi) );
199 h_assocphi_Quality075.push_back( ibooker.book1D("num_assoc(simToReco)_phi_Q075","N of associated tracks (simToReco) vs phi (Quality>0.75)",nintPhi,minPhi,maxPhi) );
200  }
201 
202  if(useLogPt){
203 BinLogX(dzres_vs_pt[j]->getTH2F());
204 BinLogX(dxyres_vs_pt[j]->getTH2F());
205 BinLogX(phires_vs_pt[j]->getTH2F());
206 BinLogX(cotThetares_vs_pt[j]->getTH2F());
207 BinLogX(ptres_vs_pt[j]->getTH2F());
208 BinLogX(h_recopT[j]->getTH1F());
209 BinLogX(h_assocpT[j]->getTH1F());
210 BinLogX(h_assoc2pT[j]->getTH1F());
211 BinLogX(h_simulpT[j]->getTH1F());
212 if (MABH) {
213 BinLogX(h_assocpT_Quality05[j]->getTH1F());
214 BinLogX(h_assocpT_Quality075[j]->getTH1F());
215 }
216  j++;
217  }
218 
219  }
220  }
221 }
222 
224  using namespace reco;
225 
226  edm::LogInfo("MuonTrackValidator") << "\n====================================================" << "\n"
227 << "Analyzing new event" << "\n"
228 << "====================================================\n" << "\n";
229  edm::ESHandle<ParametersDefinerForTP> parametersDefinerTP;
230  setup.get<TrackAssociatorRecord>().get(parametersDefiner,parametersDefinerTP);
231 
232  edm::Handle<TrackingParticleCollection> TPCollectionHeff ;
233  event.getByToken(tp_effic_Token,TPCollectionHeff);
234  const TrackingParticleCollection tPCeff = *(TPCollectionHeff.product());
235 
236  edm::Handle<TrackingParticleCollection> TPCollectionHfake ;
237  event.getByToken(tp_fake_Token,TPCollectionHfake);
238  const TrackingParticleCollection tPCfake = *(TPCollectionHfake.product());
239 
240  edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
241  event.getByToken(bsSrc_Token,recoBeamSpotHandle);
242  reco::BeamSpot bs = *recoBeamSpotHandle;
243 
244  std::vector<const reco::TrackToTrackingParticleAssociator*> associator;
245  if (UseAssociators) {
247  for (unsigned int w=0;w<associators.size();w++) {
248  event.getByLabel(associators[w],theAssociator);
249  associator.push_back( theAssociator.product() );
250  }
251  }
252 
253  int w=0;
254  for (unsigned int ww=0;ww<associators.size();ww++){
255  for (unsigned int www=0;www<label.size();www++){
256  //
257  //get collections from the event
258  //
260 
261  reco::RecoToSimCollection recSimColl;
262  reco::SimToRecoCollection simRecColl;
263  unsigned int trackCollectionSize = 0;
264 
265  if(!event.getByToken(track_Collection_Token[www], trackCollection)&&ignoremissingtkcollection_) {
266 
267 recSimColl.post_insert();
268 simRecColl.post_insert();
269 
270  }
271 
272  else {
273 
274 trackCollectionSize = trackCollection->size();
275 //associate tracks
276 if(UseAssociators){
277 edm::LogVerbatim("MuonTrackValidator") << "Analyzing "
278 << label[www].process()<<":"
279 << label[www].label()<<":"
280 << label[www].instance()<<" with "
281 << associators[ww].c_str() <<"\n";
282 
283 LogTrace("MuonTrackValidator") << "Calling associateRecoToSim method" << "\n";
284 recSimColl=associator[ww]->associateRecoToSim(trackCollection,
285  TPCollectionHfake);
286 LogTrace("MuonTrackValidator") << "Calling associateSimToReco method" << "\n";
287 simRecColl=associator[ww]->associateSimToReco(trackCollection,
288  TPCollectionHeff);
289 }
290 else{
291 edm::LogVerbatim("MuonTrackValidator") << "Analyzing "
292 << label[www].process()<<":"
293 << label[www].label()<<":"
294 << label[www].instance()<<" with "
295 << associatormap.process()<<":"
296 << associatormap.label()<<":"
297 << associatormap.instance()<<"\n";
298 
299 Handle<reco::SimToRecoCollection > simtorecoCollectionH;
300 event.getByToken(simToRecoCollection_Token,simtorecoCollectionH);
301 simRecColl= *(simtorecoCollectionH.product());
302 
303 Handle<reco::RecoToSimCollection > recotosimCollectionH;
304 event.getByToken(recoToSimCollection_Token,recotosimCollectionH);
305 recSimColl= *(recotosimCollectionH.product());
306 }
307 
308  }
309 
310 
311  //
312  //fill simulation histograms
313  //compute number of tracks per eta interval
314  //
315  edm::LogVerbatim("MuonTrackValidator") << "\n# of TrackingParticles: " << tPCeff.size() << "\n";
316  int ats = 0;
317  int st=0;
318  for (TrackingParticleCollection::size_type i=0; i<tPCeff.size(); i++){
319 bool TP_is_matched = false;
320 double quality = 0.;
321 bool Quality05 = false;
322 bool Quality075 = false;
323 
324 TrackingParticleRef tpr(TPCollectionHeff, i);
325 TrackingParticle* tp=const_cast<TrackingParticle*>(tpr.get());
326 TrackingParticle::Vector momentumTP;
327 TrackingParticle::Point vertexTP;
328 double dxySim = 0;
329 double dzSim = 0;
330 
331 //If the TrackingParticle is collison like, get the momentum and vertex at production state
332 if(parametersDefiner=="LhcParametersDefinerForTP")
333 {
334 if(! tpSelector(*tp)) continue;
335 momentumTP = tp->momentum();
336 vertexTP = tp->vertex();
337 //Calcualte the impact parameters w.r.t. PCA
338 TrackingParticle::Vector momentum = parametersDefinerTP->momentum(event,setup,tpr);
339 TrackingParticle::Point vertex = parametersDefinerTP->vertex(event,setup,tpr);
340 dxySim = (-vertex.x()*sin(momentum.phi())+vertex.y()*cos(momentum.phi()));
341 dzSim = vertex.z() - (vertex.x()*momentum.x()+vertex.y()*momentum.y())/sqrt(momentum.perp2()) * momentum.z()/sqrt(momentum.perp2());
342 }
343 //If the TrackingParticle is comics, get the momentum and vertex at PCA
344 if(parametersDefiner=="CosmicParametersDefinerForTP")
345 {
346 if(! cosmictpSelector(tpr,&bs,event,setup)) continue;
347 momentumTP = parametersDefinerTP->momentum(event,setup,tpr);
348 vertexTP = parametersDefinerTP->vertex(event,setup,tpr);
349 dxySim = (-vertexTP.x()*sin(momentumTP.phi())+vertexTP.y()*cos(momentumTP.phi()));
350 dzSim = vertexTP.z() - (vertexTP.x()*momentumTP.x()+vertexTP.y()*momentumTP.y())/sqrt(momentumTP.perp2()) * momentumTP.z()/sqrt(momentumTP.perp2());
351 }
352 edm::LogVerbatim("MuonTrackValidator") <<"--------------------Selected TrackingParticle #"<<tpr.key();
353 st++;
354 
355 h_ptSIM[w]->Fill(sqrt(momentumTP.perp2()));
356 h_etaSIM[w]->Fill(momentumTP.eta());
357 h_vertposSIM[w]->Fill(sqrt(vertexTP.perp2()));
358 
359 std::vector<std::pair<RefToBase<Track>, double> > rt;
360 if(simRecColl.find(tpr) != simRecColl.end()){
361 rt = (std::vector<std::pair<RefToBase<Track>, double> >) simRecColl[tpr];
362 if (rt.size()!=0) {
363 RefToBase<Track> assoc_recoTrack = rt.begin()->first;
364 edm::LogVerbatim("MuonTrackValidator")<<"-----------------------------associated Track #"<<assoc_recoTrack.key();
365 TP_is_matched = true;
366 ats++;
367 quality = rt.begin()->second;
368 edm::LogVerbatim("MuonTrackValidator") << "TrackingParticle #" <<tpr.key()
369 << " with pt=" << sqrt(momentumTP.perp2())
370 << " associated with quality:" << quality <<"\n";
371 if (MABH) {
372 if (quality > 0.75) {
373 Quality075 = true;
374 Quality05 = true;
375 }
376 else if (quality > 0.5) {
377 Quality05 = true;
378 }
379 }
380 }
381 }else{
382 edm::LogVerbatim("MuonTrackValidator")
383 << "TrackingParticle #" << tpr.key()
384 << " with pt,eta,phi: "
385 << sqrt(momentumTP.perp2()) << " , "
386 << momentumTP.eta() << " , "
387 << momentumTP.phi() << " , "
388 << " NOT associated to any reco::Track" << "\n";
389 }
390 
391 for (unsigned int f=0; f<etaintervals[w].size()-1; f++){
392 if (getEta(momentumTP.eta())>etaintervals[w][f]&&
393 getEta(momentumTP.eta())<etaintervals[w][f+1]) {
394 totSIMeta[w][f]++;
395 if (TP_is_matched) {
396 totASSeta[w][f]++;
397 
398 if (MABH) {
399 if (Quality075) {
400 totASSeta_Quality075[w][f]++;
401 totASSeta_Quality05[w][f]++;
402 }
403 else if (Quality05) {
404 totASSeta_Quality05[w][f]++;
405 }
406 }
407 }
408 }
409 } // END for (unsigned int f=0; f<etaintervals[w].size()-1; f++){
410 
411 for (unsigned int f=0; f<phiintervals[w].size()-1; f++){
412 if (momentumTP.phi() > phiintervals[w][f]&&
413 momentumTP.phi() <phiintervals[w][f+1]) {
414 totSIM_phi[w][f]++;
415 if (TP_is_matched) {
416 totASS_phi[w][f]++;
417 
418 if (MABH) {
419 if (Quality075) {
420 totASS_phi_Quality075[w][f]++;
421 totASS_phi_Quality05[w][f]++;
422 }
423 else if (Quality05) {
424 totASS_phi_Quality05[w][f]++;
425 }
426 }
427 }
428 }
429 } // END for (unsigned int f=0; f<phiintervals[w].size()-1; f++){
430 
431 
432 for (unsigned int f=0; f<pTintervals[w].size()-1; f++){
433  if (getPt(sqrt(momentumTP.perp2()))>pTintervals[w][f]&&
434  getPt(sqrt(momentumTP.perp2()))<pTintervals[w][f+1]) {
435  totSIMpT[w][f]++;
436 if (TP_is_matched) {
437 totASSpT[w][f]++;
438 
439 if (MABH) {
440 if (Quality075) {
441 totASSpT_Quality075[w][f]++;
442 totASSpT_Quality05[w][f]++;
443 }
444 else if (Quality05) {
445 totASSpT_Quality05[w][f]++;
446 }
447 }
448 }
449 }
450 } // END for (unsigned int f=0; f<pTintervals[w].size()-1; f++){
451 
452  for (unsigned int f=0; f<dxyintervals[w].size()-1; f++){
453 if (dxySim>dxyintervals[w][f]&&
454 dxySim<dxyintervals[w][f+1]) {
455 totSIM_dxy[w][f]++;
456 if (TP_is_matched) {
457 totASS_dxy[w][f]++;
458 }
459 }
460 } // END for (unsigned int f=0; f<dxyintervals[w].size()-1; f++){
461 
462 for (unsigned int f=0; f<dzintervals[w].size()-1; f++){
463 if (dzSim>dzintervals[w][f]&&
464 dzSim<dzintervals[w][f+1]) {
465 totSIM_dz[w][f]++;
466 if (TP_is_matched) {
467  totASS_dz[w][f]++;
468  }
469  }
470  } // END for (unsigned int f=0; f<dzintervals[w].size()-1; f++){
471 
472 for (unsigned int f=0; f<vertposintervals[w].size()-1; f++){
473 if (sqrt(vertexTP.perp2())>vertposintervals[w][f]&&
474 sqrt(vertexTP.perp2())<vertposintervals[w][f+1]) {
475 totSIM_vertpos[w][f]++;
476 if (TP_is_matched) {
477 totASS_vertpos[w][f]++;
478 }
479 }
480 } // END for (unsigned int f=0; f<vertposintervals[w].size()-1; f++){
481 
482 for (unsigned int f=0; f<zposintervals[w].size()-1; f++){
483 if (vertexTP.z()>zposintervals[w][f]&&
484 vertexTP.z()<zposintervals[w][f+1]) {
485 totSIM_zpos[w][f]++;
486 if (TP_is_matched) {
487  totASS_zpos[w][f]++;
488  }
489  }
490  } // END for (unsigned int f=0; f<zposintervals[w].size()-1; f++){
491 
492 int nSimHits = 0;
493 if (usetracker && usemuon) {
494 nSimHits= tpr.get()->numberOfHits();
495 }
496 else if (!usetracker && usemuon) {
497 nSimHits= tpr.get()->numberOfHits() - tpr.get()->numberOfTrackerHits();
498 }
499 else if (usetracker && !usemuon) {
500 nSimHits=tpr.get()->numberOfTrackerHits();
501 }
502 
503 
504  int tmp = std::min(nSimHits,int(maxHit-1));
505 edm::LogVerbatim("MuonTrackValidator") << "\t N simhits = "<< nSimHits<<"\n";
506 
507 totSIM_hit[w][tmp]++;
508 if (TP_is_matched) totASS_hit[w][tmp]++;
509 
510 if (TP_is_matched)
511 {
512 RefToBase<Track> assoctrack = rt.begin()->first;
513 nrecHit_vs_nsimHit_sim2rec[w]->Fill( assoctrack->numberOfValidHits(),nSimHits);
514 }
515  } // End for (TrackingParticleCollection::size_type i=0; i<tPCeff.size(); i++){
516  if (st!=0) h_tracksSIM[w]->Fill(st);
517 
518 
519  //
520  //fill reconstructed track histograms
521  //
522  edm::LogVerbatim("MuonTrackValidator") << "\n# of reco::Tracks with "
523 << label[www].process()<<":"
524 << label[www].label()<<":"
525 << label[www].instance()
526 << ": " << trackCollectionSize << "\n";
527  int at=0;
528  int rT=0;
529  for(edm::View<Track>::size_type i=0; i<trackCollectionSize; ++i){
530  bool Track_is_matched = false;
531 RefToBase<Track> track(trackCollection, i);
532 rT++;
533 
534 std::vector<std::pair<TrackingParticleRef, double> > tp;
536 
537 // new logic (bidirectional)
538 if (BiDirectional_RecoToSim_association) {
539 edm::LogVerbatim("MuonTrackValidator")<<"----------------------------------------Track #"<< track.key();
540 
541 if(recSimColl.find(track) != recSimColl.end()) {
542 tp = recSimColl[track];
543 if (tp.size() != 0) {
544 tpr = tp.begin()->first;
545 // RtS and StR must associate the same pair !
546 if(simRecColl.find(tpr) != simRecColl.end()) {
547 std::vector<std::pair<RefToBase<Track>, double> > track_checkback = simRecColl[tpr];
548 RefToBase<Track> assoc_track_checkback;
549 assoc_track_checkback = track_checkback.begin()->first;
550 
551 if ( assoc_track_checkback.key() == track.key() ) {
552 edm::LogVerbatim("MuonTrackValidator")<<"------------------associated TrackingParticle #"<<tpr.key();
553 Track_is_matched = true;
554 at++;
555 double Purity = tp.begin()->second;
556 double Quality = track_checkback.begin()->second;
557 edm::LogVerbatim("MuonTrackValidator") << "reco::Track #" << track.key() << " with pt=" << track->pt()
558 << " associated with quality:" << Purity <<"\n";
559 if (MABH) h_PurityVsQuality[w]->Fill(Quality,Purity);
560 }
561 }
562 }
563 }
564 
565 if (!Track_is_matched) edm::LogVerbatim("MuonTrackValidator")
566 << "reco::Track #" << track.key() << " with pt=" << track->pt() << " NOT associated to any TrackingParticle" << "\n";
567 }
568 // old logic (bugged)
569 else {
570 if(recSimColl.find(track) != recSimColl.end()){
571 tp = recSimColl[track];
572 if (tp.size()!=0) {
573 Track_is_matched = true;
574 tpr = tp.begin()->first;
575 at++;
576 edm::LogVerbatim("MuonTrackValidator") << "reco::Track #" << track.key() << " with pt=" << track->pt()
577 << " associated with quality:" << tp.begin()->second <<"\n";
578 }
579 } else {
580 edm::LogVerbatim("MuonTrackValidator") << "reco::Track #" << track.key() << " with pt=" << track->pt()
581 << " NOT associated to any TrackingParticle" << "\n";
582 }
583 }
584 
585 //Compute fake rate vs eta
586 for (unsigned int f=0; f<etaintervals[w].size()-1; f++){
587 if (getEta(track->momentum().eta())>etaintervals[w][f]&&
588 getEta(track->momentum().eta())<etaintervals[w][f+1]) {
589 totRECeta[w][f]++;
590 if (Track_is_matched) {
591 totASS2eta[w][f]++;
592 }
593 }
594 } // End for (unsigned int f=0; f<etaintervals[w].size()-1; f++){
595 
596 for (unsigned int f=0; f<phiintervals[w].size()-1; f++){
597 if (track->momentum().phi()>phiintervals[w][f]&&
598 track->momentum().phi()<phiintervals[w][f+1]) {
599 totREC_phi[w][f]++;
600 if (Track_is_matched) {
601 totASS2_phi[w][f]++;
602 }
603 }
604 } // End for (unsigned int f=0; f<phiintervals[w].size()-1; f++){
605 
606 
607 for (unsigned int f=0; f<pTintervals[w].size()-1; f++){
608 if (getPt(sqrt(track->momentum().perp2()))>pTintervals[w][f]&&
609 getPt(sqrt(track->momentum().perp2()))<pTintervals[w][f+1]) {
610 totRECpT[w][f]++;
611 if (Track_is_matched) {
612 totASS2pT[w][f]++;
613 }
614 }
615 } // End for (unsigned int f=0; f<pTintervals[w].size()-1; f++){
616 
617 for (unsigned int f=0; f<dxyintervals[w].size()-1; f++){
618 if (track->dxy(bs.position())>dxyintervals[w][f]&&
619 track->dxy(bs.position())<dxyintervals[w][f+1]) {
620 totREC_dxy[w][f]++;
621 if (Track_is_matched) {
622 totASS2_dxy[w][f]++;
623 }
624 }
625 } // End for (unsigned int f=0; f<dxyintervals[w].size()-1; f++){
626 
627 for (unsigned int f=0; f<dzintervals[w].size()-1; f++){
628 if (track->dz(bs.position())>dzintervals[w][f]&&
629 track->dz(bs.position())<dzintervals[w][f+1]) {
630 totREC_dz[w][f]++;
631 if (Track_is_matched) {
632 totASS2_dz[w][f]++;
633 }
634 }
635 } // End for (unsigned int f=0; f<dzintervals[w].size()-1; f++){
636 
637 int tmp = std::min((int)track->found(),int(maxHit-1));
638  totREC_hit[w][tmp]++;
639 if (Track_is_matched) totASS2_hit[w][tmp]++;
640 
641 edm::LogVerbatim("MuonTrackValidator") << "\t N valid rechits = "<< (int)track->found() <<"\n";
642 
643 //Fill other histos
644  try{
645 if (!Track_is_matched) continue;
646 
647 if (associators[ww]=="trackAssociatorByChi2"){
648 //association chi2
649 double assocChi2 = -tp.begin()->second;//in association map is stored -chi2
650 h_assochi2[www]->Fill(assocChi2);
651 h_assochi2_prob[www]->Fill(TMath::Prob((assocChi2)*5,5));
652 }
653 else if (associators[ww]=="trackAssociatorByHits"){
654 double fraction = tp.begin()->second;
655 h_assocFraction[www]->Fill(fraction);
656 h_assocSharedHit[www]->Fill(fraction*track->numberOfValidHits());
657 }
658 
659 //nchi2 and hits global distributions
660 h_nchi2[w]->Fill(track->normalizedChi2());
661 h_nchi2_prob[w]->Fill(TMath::Prob(track->chi2(),(int)track->ndof()));
662 h_hits[w]->Fill(track->numberOfValidHits());
663 h_losthits[w]->Fill(track->numberOfLostHits());
664 chi2_vs_nhits[w]->Fill(track->numberOfValidHits(),track->normalizedChi2());
665 h_charge[w]->Fill( track->charge() );
666 
667 //Get tracking particle parameters at point of closest approach to the beamline
668 TrackingParticle::Vector momentumTP = parametersDefinerTP->momentum(event,setup,tpr) ;
669 TrackingParticle::Point vertexTP = parametersDefinerTP->vertex(event,setup,tpr);
670 double ptSim = sqrt(momentumTP.perp2());
671 double qoverpSim = tpr->charge()/sqrt(momentumTP.x()*momentumTP.x()+momentumTP.y()*momentumTP.y()+momentumTP.z()*momentumTP.z());
672 double thetaSim = momentumTP.theta();
673 double lambdaSim = M_PI/2-momentumTP.theta();
674 double phiSim = momentumTP.phi();
675 double dxySim = (-vertexTP.x()*sin(momentumTP.phi())+vertexTP.y()*cos(momentumTP.phi()));
676 double dzSim = vertexTP.z() - (vertexTP.x()*momentumTP.x()+vertexTP.y()*momentumTP.y())/sqrt(momentumTP.perp2()) * momentumTP.z()/sqrt(momentumTP.perp2());
677 
678 // removed unused variable, left this in case it has side effects
679 track->parameters();
680 
681 double qoverpRec(0);
682 double qoverpErrorRec(0);
683 double ptRec(0);
684 double ptErrorRec(0);
685 double lambdaRec(0);
686 double lambdaErrorRec(0);
687 double phiRec(0);
688 double phiErrorRec(0);
689 
690 
691 //loop to decide whether to take gsfTrack (utilisation of mode-function) or common track
692 const GsfTrack* gsfTrack(0);
693 if(useGsf){
694 gsfTrack = dynamic_cast<const GsfTrack*>(&(*track));
695 if (gsfTrack==0) edm::LogInfo("MuonTrackValidator") << "Trying to access mode for a non-GsfTrack";
696 }
697 
698 if (gsfTrack) {
699 // get values from mode
700 getRecoMomentum(*gsfTrack, ptRec, ptErrorRec, qoverpRec, qoverpErrorRec,
701 lambdaRec,lambdaErrorRec, phiRec, phiErrorRec);
702 }
703 
704 else {
705 // get values from track (without mode)
706 getRecoMomentum(*track, ptRec, ptErrorRec, qoverpRec, qoverpErrorRec,
707 lambdaRec,lambdaErrorRec, phiRec, phiErrorRec);
708 }
709 
710 double thetaRec = track->theta();
711 double ptError = ptErrorRec;
712 double ptres = ptRec - ptSim;
713 double etares = track->eta()-momentumTP.Eta();
714 double dxyRec = track->dxy(bs.position());
715 double dzRec = track->dz(bs.position());
716 // eta residue; pt, k, theta, phi, dxy, dz pulls
717 double qoverpPull=(qoverpRec-qoverpSim)/qoverpErrorRec;
718 double thetaPull=(lambdaRec-lambdaSim)/lambdaErrorRec;
719 double phiDiff = phiRec - phiSim;
720 if (abs(phiDiff) > M_PI) {
721 if (phiDiff >0.) phiDiff = phiDiff - 2.*M_PI;
722  else phiDiff = phiDiff + 2.*M_PI;
723 }
724 double phiPull=phiDiff/phiErrorRec;
725 double dxyPull=(dxyRec-dxySim)/track->dxyError();
726 double dzPull=(dzRec-dzSim)/track->dzError();
727 
728 double contrib_Qoverp = ((qoverpRec-qoverpSim)/qoverpErrorRec)*
729 ((qoverpRec-qoverpSim)/qoverpErrorRec)/5;
730 double contrib_dxy = ((dxyRec-dxySim)/track->dxyError())*((dxyRec-dxySim)/track->dxyError())/5;
731 double contrib_dz = ((dzRec-dzSim)/track->dzError())*((dzRec-dzSim)/track->dzError())/5;
732 double contrib_theta = ((lambdaRec-lambdaSim)/lambdaErrorRec)*
733 ((lambdaRec-lambdaSim)/lambdaErrorRec)/5;
734 double contrib_phi = (phiDiff/phiErrorRec)*(phiDiff/phiErrorRec)/5;
735 
736 LogTrace("MuonTrackValidator") << "assocChi2=" << tp.begin()->second << "\n"
737 << "" << "\n"
738 << "ptREC=" << ptRec << "\n"
739 << "etaREC=" << track->eta() << "\n"
740 << "qoverpREC=" << qoverpRec << "\n"
741 << "dxyREC=" << dxyRec << "\n"
742 << "dzREC=" << dzRec << "\n"
743 << "thetaREC=" << track->theta() << "\n"
744 << "phiREC=" << phiRec << "\n"
745 << "" << "\n"
746 << "qoverpError()=" << qoverpErrorRec << "\n"
747 << "dxyError()=" << track->dxyError() << "\n"
748 << "dzError()=" << track->dzError() << "\n"
749 << "thetaError()=" << lambdaErrorRec << "\n"
750 << "phiError()=" << phiErrorRec << "\n"
751 << "" << "\n"
752 << "ptSIM=" << ptSim << "\n"
753 << "etaSIM=" << momentumTP.Eta() << "\n"
754 << "qoverpSIM=" << qoverpSim << "\n"
755 << "dxySIM=" << dxySim << "\n"
756 << "dzSIM=" << dzSim << "\n"
757 << "thetaSIM=" << M_PI/2-lambdaSim << "\n"
758 << "phiSIM=" << phiSim << "\n"
759 << "" << "\n"
760 << "contrib_Qoverp=" << contrib_Qoverp << "\n"
761 << "contrib_dxy=" << contrib_dxy << "\n"
762 << "contrib_dz=" << contrib_dz << "\n"
763 << "contrib_theta=" << contrib_theta << "\n"
764 << "contrib_phi=" << contrib_phi << "\n"
765 << "" << "\n"
766 <<"chi2PULL="<<contrib_Qoverp+contrib_dxy+contrib_dz+contrib_theta+contrib_phi<<"\n";
767 
768 h_pullQoverp[w]->Fill(qoverpPull);
769 h_pullTheta[w]->Fill(thetaPull);
770 h_pullPhi[w]->Fill(phiPull);
771 h_pullDxy[w]->Fill(dxyPull);
772 h_pullDz[w]->Fill(dzPull);
773 
774 
775 h_pt[w]->Fill(ptres/ptError);
776 h_eta[w]->Fill(etares);
777 etares_vs_eta[w]->Fill(getEta(track->eta()),etares);
778 
779 
780 //chi2 and #hit vs eta: fill 2D histos
781 chi2_vs_eta[w]->Fill(getEta(track->eta()),track->normalizedChi2());
782 nhits_vs_eta[w]->Fill(getEta(track->eta()),track->numberOfValidHits());
783 nDThits_vs_eta[w]->Fill(getEta(track->eta()),track->hitPattern().numberOfValidMuonDTHits());
784 nCSChits_vs_eta[w]->Fill(getEta(track->eta()),track->hitPattern().numberOfValidMuonCSCHits());
785 nRPChits_vs_eta[w]->Fill(getEta(track->eta()),track->hitPattern().numberOfValidMuonRPCHits());
786 
787 nlosthits_vs_eta[w]->Fill(getEta(track->eta()),track->numberOfLostHits());
788 
789 //resolution of track params: fill 2D histos
790 dxyres_vs_eta[w]->Fill(getEta(track->eta()),dxyRec-dxySim);
791 ptres_vs_eta[w]->Fill(getEta(track->eta()),(ptRec-ptSim)/ptRec);
792 dzres_vs_eta[w]->Fill(getEta(track->eta()),dzRec-dzSim);
793 phires_vs_eta[w]->Fill(getEta(track->eta()),phiDiff);
794 cotThetares_vs_eta[w]->Fill(getEta(track->eta()), cos(thetaRec)/sin(thetaRec) - cos(thetaSim)/sin(thetaSim));
795 
796 //same as before but vs pT
797 dxyres_vs_pt[w]->Fill(getPt(ptRec),dxyRec-dxySim);
798 ptres_vs_pt[w]->Fill(getPt(ptRec),(ptRec-ptSim)/ptRec);
799 dzres_vs_pt[w]->Fill(getPt(ptRec),dzRec-dzSim);
800 phires_vs_pt[w]->Fill(getPt(ptRec),phiDiff);
801  cotThetares_vs_pt[w]->Fill(getPt(ptRec), cos(thetaRec)/sin(thetaRec) - cos(thetaSim)/sin(thetaSim));
802 
803 //pulls of track params vs eta: fill 2D histos
804 dxypull_vs_eta[w]->Fill(getEta(track->eta()),dxyPull);
805 ptpull_vs_eta[w]->Fill(getEta(track->eta()),ptres/ptError);
806 dzpull_vs_eta[w]->Fill(getEta(track->eta()),dzPull);
807 phipull_vs_eta[w]->Fill(getEta(track->eta()),phiPull);
808 thetapull_vs_eta[w]->Fill(getEta(track->eta()),thetaPull);
809 
810 //plots vs phi
811 nhits_vs_phi[w]->Fill(phiRec,track->numberOfValidHits());
812 chi2_vs_phi[w]->Fill(phiRec,track->normalizedChi2());
813 ptmean_vs_eta_phi[w]->Fill(phiRec,getEta(track->eta()),ptRec);
814 phimean_vs_eta_phi[w]->Fill(phiRec,getEta(track->eta()),phiRec);
815 ptres_vs_phi[w]->Fill(phiRec,(ptRec-ptSim)/ptRec);
816 phires_vs_phi[w]->Fill(phiRec,phiDiff);
817 ptpull_vs_phi[w]->Fill(phiRec,ptres/ptError);
818 phipull_vs_phi[w]->Fill(phiRec,phiPull);
819 thetapull_vs_phi[w]->Fill(phiRec,thetaPull);
820 
821 int nSimHits = 0;
822 if (usetracker && usemuon) {
823 nSimHits= tpr.get()->numberOfHits();
824 }
825 else if (!usetracker && usemuon) {
826 nSimHits= tpr.get()->numberOfHits() - tpr.get()->numberOfTrackerHits();
827 }
828 else if (usetracker && !usemuon) {
829 nSimHits=tpr.get()->numberOfTrackerHits();
830 }
831 
832 nrecHit_vs_nsimHit_rec2sim[w]->Fill(track->numberOfValidHits(), nSimHits);
833 
834 } // End of try{
835 catch (cms::Exception e){
836 LogTrace("MuonTrackValidator") << "exception found: " << e.what() << "\n";
837 }
838  } // End of for(View<Track>::size_type i=0; i<trackCollectionSize; ++i){
839  if (at!=0) h_tracks[w]->Fill(at);
840  h_fakes[w]->Fill(rT-at);
841  edm::LogVerbatim("MuonTrackValidator") << "Total Simulated: " << st << "\n"
842 << "Total Associated (simToReco): " << ats << "\n"
843 << "Total Reconstructed: " << rT << "\n"
844 << "Total Associated (recoToSim): " << at << "\n"
845 << "Total Fakes: " << rT-at << "\n";
846  nrec_vs_nsim[w]->Fill(rT,st);
847  w++;
848  } // End of for (unsigned int www=0;www<label.size();www++){
849  } //END of for (unsigned int ww=0;ww<associators.size();ww++){
850 }
851 
853 
854  int w=0;
855  for (unsigned int ww=0;ww<associators.size();ww++){
856  for (unsigned int www=0;www<label.size();www++){
857 
858  //chi2 and #hit vs eta: get mean from 2D histos
859  doProfileX(chi2_vs_eta[w],h_chi2meanh[w]);
860  doProfileX(nhits_vs_eta[w],h_hits_eta[w]);
861  doProfileX(nDThits_vs_eta[w],h_DThits_eta[w]);
862  doProfileX(nCSChits_vs_eta[w],h_CSChits_eta[w]);
863  doProfileX(nRPChits_vs_eta[w],h_RPChits_eta[w]);
864 
865  doProfileX(nlosthits_vs_eta[w],h_losthits_eta[w]);
866  //vs phi
867  doProfileX(chi2_vs_nhits[w],h_chi2meanhitsh[w]);
868  doProfileX(chi2_vs_phi[w],h_chi2mean_vs_phi[w]);
869  doProfileX(nhits_vs_phi[w],h_hits_phi[w]);
870 
871  fillPlotFromVector(h_recoeta[w],totRECeta[w]);
872  fillPlotFromVector(h_simuleta[w],totSIMeta[w]);
873  fillPlotFromVector(h_assoceta[w],totASSeta[w]);
874  fillPlotFromVector(h_assoc2eta[w],totASS2eta[w]);
875 
876  fillPlotFromVector(h_recopT[w],totRECpT[w]);
877  fillPlotFromVector(h_simulpT[w],totSIMpT[w]);
878  fillPlotFromVector(h_assocpT[w],totASSpT[w]);
879  fillPlotFromVector(h_assoc2pT[w],totASS2pT[w]);
880 
881  fillPlotFromVector(h_recohit[w],totREC_hit[w]);
882  fillPlotFromVector(h_simulhit[w],totSIM_hit[w]);
883  fillPlotFromVector(h_assochit[w],totASS_hit[w]);
884  fillPlotFromVector(h_assoc2hit[w],totASS2_hit[w]);
885 
886  fillPlotFromVector(h_recophi[w],totREC_phi[w]);
887  fillPlotFromVector(h_simulphi[w],totSIM_phi[w]);
888  fillPlotFromVector(h_assocphi[w],totASS_phi[w]);
889  fillPlotFromVector(h_assoc2phi[w],totASS2_phi[w]);
890 
891  fillPlotFromVector(h_recodxy[w],totREC_dxy[w]);
892  fillPlotFromVector(h_simuldxy[w],totSIM_dxy[w]);
893  fillPlotFromVector(h_assocdxy[w],totASS_dxy[w]);
894  fillPlotFromVector(h_assoc2dxy[w],totASS2_dxy[w]);
895 
896  fillPlotFromVector(h_recodz[w],totREC_dz[w]);
897  fillPlotFromVector(h_simuldz[w],totSIM_dz[w]);
898  fillPlotFromVector(h_assocdz[w],totASS_dz[w]);
899  fillPlotFromVector(h_assoc2dz[w],totASS2_dz[w]);
900 
901  fillPlotFromVector(h_simulvertpos[w],totSIM_vertpos[w]);
902  fillPlotFromVector(h_assocvertpos[w],totASS_vertpos[w]);
903 
904  fillPlotFromVector(h_simulzpos[w],totSIM_zpos[w]);
905  fillPlotFromVector(h_assoczpos[w],totASS_zpos[w]);
906 
907  if (MABH) {
908 fillPlotFromVector(h_assoceta_Quality05[w] ,totASSeta_Quality05[w]);
909 fillPlotFromVector(h_assoceta_Quality075[w],totASSeta_Quality075[w]);
910 fillPlotFromVector(h_assocpT_Quality05[w] ,totASSpT_Quality05[w]);
911 fillPlotFromVector(h_assocpT_Quality075[w],totASSpT_Quality075[w]);
912 fillPlotFromVector(h_assocphi_Quality05[w] ,totASS_phi_Quality05[w]);
913 fillPlotFromVector(h_assocphi_Quality075[w],totASS_phi_Quality075[w]);
914  }
915 
916  w++;
917  }
918  }
919 
920  if ( out.size() != 0 && dbe_ ) dbe_->save(out);
921 }
922 
923 
924 void
925 MuonTrackValidator::getRecoMomentum (const reco::Track& track, double& pt, double& ptError,
926 double& qoverp, double& qoverpError, double& lambda,double& lambdaError, double& phi, double& phiError ) const {
927  pt = track.pt();
928  ptError = track.ptError();
929  qoverp = track.qoverp();
930  qoverpError = track.qoverpError();
931  lambda = track.lambda();
932  lambdaError = track.lambdaError();
933  phi = track.phi();
934  phiError = track.phiError();
935 
936 }
937 
938 void
939 MuonTrackValidator::getRecoMomentum (const reco::GsfTrack& gsfTrack, double& pt, double& ptError,
940 double& qoverp, double& qoverpError, double& lambda,double& lambdaError, double& phi, double& phiError ) const {
941 
942  pt = gsfTrack.ptMode();
943  ptError = gsfTrack.ptModeError();
944  qoverp = gsfTrack.qoverpMode();
945  qoverpError = gsfTrack.qoverpModeError();
946  lambda = gsfTrack.lambdaMode();
947  lambdaError = gsfTrack.lambdaModeError();
948  phi = gsfTrack.phiMode();
949  phiError = gsfTrack.phiModeError();
950 
951 }
double qoverp() const
q / p
Definition: TrackBase.h:621
virtual char const * what() const
Definition: Exception.cc:141
double phiModeError() const
error on phi from mode
Definition: GsfTrack.h:94
unsigned int size_type
Definition: View.h:85
int i
Definition: DBlmapReader.cc:9
const double w
Definition: UKUtility.cc:23
std::vector< TrackingParticle > TrackingParticleCollection
double lambdaMode() const
Lambda angle from mode.
Definition: GsfTrack.h:45
const_iterator end() const
last iterator over the map (read only)
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:609
MonitorElement * bookProfile(Args &&...args)
Definition: DQMStore.h:157
tuple maxDxy
Definition: gather_cfg.py:52
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:449
void cd(void)
Definition: DQMStore.cc:266
double theta() const
polar angle
Definition: TrackBase.h:627
double dxyError() const
error on dxy
Definition: TrackBase.h:844
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
const_iterator find(const key_type &k) const
find element with specified reference key
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:693
unsigned short numberOfLostHits() const
number of cases where track crossed a layer without getting a hit.
Definition: TrackBase.h:874
key_type key() const
Accessor for product key.
Definition: Ref.h:266
const Vector & momentum() const
track momentum vector
Definition: TrackBase.h:723
void getRecoMomentum(const reco::Track &track, double &pt, double &ptError, double &qoverp, double &qoverpError, double &lambda, double &lambdaError, double &phi, double &phiError) const
retrieval of reconstructed momentum components from reco::Track (== mean values for GSF) ...
uint16_t size_type
int numberOfValidMuonCSCHits() const
Definition: HitPattern.h:794
double qoverpMode() const
q/p from mode
Definition: GsfTrack.h:41
int numberOfValidMuonRPCHits() const
Definition: HitPattern.h:799
double ptModeError() const
error on Pt (set to 1000 TeV if charge==0 for safety) from mode
Definition: GsfTrack.h:81
math::XYZPointD Point
point in the space
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:699
double chi2() const
chi-squared of the fit
Definition: TrackBase.h:597
MonitorElement * bookProfile2D(Args &&...args)
Definition: DQMStore.h:163
size_t key() const
Definition: RefToBase.h:228
void post_insert()
post insert action
double ndof() const
number of degrees of freedom of the fit
Definition: TrackBase.h:603
T sqrt(T t)
Definition: SSEVec.h:48
double pt() const
track transverse momentum
Definition: TrackBase.h:669
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:811
double phiError() const
error on phi
Definition: TrackBase.h:838
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:115
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int j
Definition: DBlmapReader.cc:9
double lambda() const
Lambda angle.
Definition: TrackBase.h:633
double f[11][100]
TH1F * getTH1F(std::string name, std::string process, std::string rootfolder, DQMStore *dbe_, bool verb, bool clone)
void analyze(const edm::Event &, const edm::EventSetup &)
Method called before the event loop.
unsigned short numberOfValidHits() const
number of valid hits found
Definition: TrackBase.h:868
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:242
T min(T a, T b)
Definition: MathUtil.h:58
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
#define LogTrace(id)
#define M_PI
DQMStore * dbe_
int nint(float a)
Return the nearest integer - analogous to the FORTRAN intrinsic NINT.
Definition: nint.h:8
double qoverpError() const
error on signed transverse curvature
Definition: TrackBase.h:805
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:657
double dzError() const
error on dz
Definition: TrackBase.h:862
tuple out
Definition: dbtoconf.py:99
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:274
void goUp(void)
Definition: DQMStore.cc:278
T const * product() const
Definition: Handle.h:81
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:133
const HitPattern & hitPattern() const
Access the hit pattern, indicating in which Tracker layers the track has hits.
Definition: TrackBase.h:384
double qoverpModeError() const
error on signed transverse curvature from mode
Definition: GsfTrack.h:79
const T & get() const
Definition: EventSetup.h:55
std::string const & label() const
Definition: InputTag.h:42
Point vertex() const
Parent vertex position.
ParameterVector parameters() const
Track parameters with one-to-one correspondence to the covariance matrix.
Definition: TrackBase.h:773
std::string const & process() const
Definition: InputTag.h:46
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
double lambdaError() const
error on lambda
Definition: TrackBase.h:826
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &)
unsigned short found() const
Number of valid hits on track.
Definition: Track.h:194
int numberOfValidMuonDTHits() const
Definition: HitPattern.h:789
Vector momentum() const
spatial momentum vector
TH2F * getTH2F(std::string name, std::string process, std::string rootfolder, DQMStore *dbe_, bool verb, bool clone)
double phiMode() const
azimuthal angle of momentum vector from mode
Definition: GsfTrack.h:57
Monte Carlo truth information used for tracking validation.
int charge() const
track electric charge
Definition: TrackBase.h:615
double lambdaModeError() const
error on lambda from mode
Definition: GsfTrack.h:90
math::XYZVectorD Vector
point in the space
double dxy() const
dxy parameter. (This is the transverse impact parameter w.r.t. to (0,0,0) ONLY if refPoint is close t...
Definition: TrackBase.h:639
double ptMode() const
track transverse momentum from mode
Definition: GsfTrack.h:49
std::string const & instance() const
Definition: InputTag.h:43
void setup(std::vector< TH2F > &depth, std::string name, std::string units="")
Definition: Run.h:41
list at
Definition: asciidump.py:428
void endRun(edm::Run const &, edm::EventSetup const &)
Method called at the end of the event loop.
Definition: DDAxes.h:10