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