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