CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TestSmoothHits.cc
Go to the documentation of this file.
2 
11 #include <TDirectory.h>
13 
21 
24 using namespace std;
25 using namespace edm;
26 
28  conf_(iConfig){
29  LogTrace("TestSmoothHits") << conf_<< std::endl;
30  propagatorName = conf_.getParameter<std::string>("Propagator");
31  builderName = conf_.getParameter<std::string>("TTRHBuilder");
32  srcName = conf_.getParameter<std::string>("src");
33  fname = conf_.getParameter<std::string>("Fitter");
34  sname = conf_.getParameter<std::string>("Smoother");
35  mineta = conf_.getParameter<double>("mineta");
36  maxeta = conf_.getParameter<double>("maxeta");
37 }
38 
40 
42 {
43 
44  iSetup.get<TrackerDigiGeometryRecord>().get(theG);
45  iSetup.get<IdealMagneticFieldRecord>().get(theMF);
48  iSetup.get<TrajectoryFitter::Record>().get(fname, fit);
49  iSetup.get<TrajectoryFitter::Record>().get(sname, smooth);
50 
51  file = new TFile("testSmoothHits.root","recreate");
52  for (int i=0; i!=6; i++)
53  for (int j=0; j!=9; j++){
54  if (i==0 && j>2) break;
55  if (i==1 && j>1) break;
56  if (i==2 && j>3) break;
57  if (i==3 && j>2) break;
58  if (i==4 && j>5) break;
59  if (i==5 && j>8) break;
60  title.str("");
61  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
62  hPullGP_X_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
63  title.str("");
64  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts";
65  hPullGP_Y_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
66  title.str("");
67  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts";
68  hPullGP_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
69  title.str("");
70  title << "Chi2Increment_" << i+1 << "-" << j+1;
71  hChi2Increment[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
72 
73  title.str("");
74  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
75  hPullGM_X_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
76  title.str("");
77  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
78  hPullGM_Y_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
79  title.str("");
80  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
81  hPullGM_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
82 
83  title.str("");
84  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
85  hPullGP_X_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
86  title.str("");
87  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
88  hPullGP_Y_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
89  title.str("");
90  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
91  hPullGP_Z_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
92 
93  title.str("");
94  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
95  hPullGP_X_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
96  title.str("");
97  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
98  hPullGP_Y_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
99  title.str("");
100  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
101  hPullGP_Z_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
102 
103  if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
104  //mono
105  title.str("");
106  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
107  hPullGP_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
108  title.str("");
109  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
110  hPullGP_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
111  title.str("");
112  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
113  hPullGP_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
114 
115  title.str("");
116  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
117  hPullGM_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
118  title.str("");
119  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
120  hPullGM_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
121  title.str("");
122  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
123  hPullGM_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
124 
125  title.str("");
126  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
127  hPullGP_X_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
128  title.str("");
129  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
130  hPullGP_Y_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
131  title.str("");
132  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
133  hPullGP_Z_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
134 
135  title.str("");
136  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
137  hPullGP_X_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
138  title.str("");
139  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
140  hPullGP_Y_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
141  title.str("");
142  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
143  hPullGP_Z_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
144 
145  //stereo
146  title.str("");
147  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
148  hPullGP_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
149  title.str("");
150  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
151  hPullGP_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
152  title.str("");
153  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
154  hPullGP_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
155 
156  title.str("");
157  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
158  hPullGM_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
159  title.str("");
160  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
161  hPullGM_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
162  title.str("");
163  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
164  hPullGM_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
165 
166  title.str("");
167  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
168  hPullGP_X_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
169  title.str("");
170  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
171  hPullGP_Y_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
172  title.str("");
173  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
174  hPullGP_Z_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
175 
176  title.str("");
177  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
178  hPullGP_X_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
179  title.str("");
180  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
181  hPullGP_Y_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
182  title.str("");
183  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
184  hPullGP_Z_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
185  }
186  }
187  hTotChi2Increment = new TH1F("TotChi2Increment","TotChi2Increment",1000,0,100);
188  hChi2_vs_Process = new TH2F("Chi2_vs_Process","Chi2_vs_Process",1000,0,100,17,-0.5,16.5);
189  hChi2_vs_clsize = new TH2F("Chi2_vs_clsize","Chi2_vs_clsize",1000,0,100,17,-0.5,16.5);
190 }
191 
192 
194 {
195  LogTrace("TestSmoothHits") << "new event" << std::endl;
196 
198  hitAssociator = new TrackerHitAssociator(iEvent);
199 
201 
202  for (TrackCandidateCollection::const_iterator i=theTCCollection->begin(); i!=theTCCollection->end();i++){
203 
204  LogTrace("TestSmoothHits") << "new candidate" << std::endl;
205 
206  const TrackCandidate * theTC = &(*i);
208  const TrackCandidate::range& recHitVec=theTC->recHits();
209 
210  //convert PTrajectoryStateOnDet to TrajectoryStateOnSurface
211 
212 
213  DetId detId(state.detId());
214  TrajectoryStateOnSurface theTSOS=
215  trajectoryStateTransform::transientState(state, &(theG->idToDet(detId)->surface()),theMF.product());
216 
217  if (theTSOS.globalMomentum().eta()>maxeta || theTSOS.globalMomentum().eta()<mineta) continue;
218 
219  //convert the TrackingRecHit vector to a TransientTrackingRecHit vector
221 
223  i!=recHitVec.second; i++){
224  hits.push_back(theBuilder->build(&(*i) ));
225  }
226 
227  //call the fitter
228  std::vector<Trajectory> fitted = fit->fit(theTC->seed(), hits, theTSOS);
229  //call the smoother
230  std::vector<Trajectory> result;
231  for(std::vector<Trajectory>::iterator it = fitted.begin(); it != fitted.end(); it++) {
232  std::vector<Trajectory> smoothed = smooth->trajectories(*it);
233  result.insert(result.end(), smoothed.begin(), smoothed.end());
234  }
235  if (result.size()==0) continue;
236  std::vector<TrajectoryMeasurement> vtm = result[0].measurements();
237 
238  TSOS lastState = theTSOS;
239  for (std::vector<TrajectoryMeasurement>::iterator tm=vtm.begin(); tm!=vtm.end();tm++){
240 
242  if ((rhit)->isValid()==0&&rhit->det()!=0) continue;
243  LogTrace("TestSmoothHits") << "new hit" ;
244 
245  int subdetId = rhit->det()->geographicalId().subdetId();
246  int layerId = 0;
247  DetId id = rhit->det()->geographicalId();
248  if (id.subdetId()==3) layerId = ((TIBDetId)(id)).layer();
249  if (id.subdetId()==5) layerId = ((TOBDetId)(id)).layer();
250  if (id.subdetId()==1) layerId = ((PXBDetId)(id)).layer();
251  if (id.subdetId()==4) layerId = ((TIDDetId)(id)).wheel();
252  if (id.subdetId()==6) layerId = ((TECDetId)(id)).wheel();
253  if (id.subdetId()==2) layerId = ((PXFDetId)(id)).disk();
254  LogTrace("TestSmoothHits") << "subdetId=" << subdetId << " layerId=" << layerId ;
255 
256  double delta = 99999;
257  LocalPoint rhitLPv = rhit->localPosition();
258 
259  std::vector<PSimHit> assSimHits = hitAssociator->associateHit(*(rhit->hit()));
260  if (assSimHits.size()==0) continue;
261  PSimHit shit;
262  for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
263  if ((m->localPosition()-rhitLPv).mag()<delta) {
264  shit=*m;
265  delta = (m->localPosition()-rhitLPv).mag();
266  }
267  }
268 
269  TSOS currentState = combiner(tm->backwardPredictedState(), tm->forwardPredictedState());
270  TSOS updatedState = tm->updatedState();
271 
272  //plot chi2 increment
273  double chi2increment = tm->estimate();
274  LogTrace("TestSmoothHits") << "tm->estimate()=" << tm->estimate();
275  title.str("");
276  title << "Chi2Increment_" << subdetId << "-" << layerId;
277  hChi2Increment[title.str()]->Fill( chi2increment );
278  hTotChi2Increment->Fill( chi2increment );
279  hChi2_vs_Process->Fill( chi2increment, shit.processType() );
280  if (dynamic_cast<const SiPixelRecHit*>(rhit->hit()))
281  hChi2_vs_clsize->Fill( chi2increment, ((const SiPixelRecHit*)(rhit->hit()))->cluster()->size() );
282  if (dynamic_cast<const SiStripRecHit2D*>(rhit->hit()))
283  hChi2_vs_clsize->Fill( chi2increment, ((const SiStripRecHit2D*)(rhit->hit()))->cluster()->amplitudes().size() );
284 
285  //test hits
286  const Surface * surf = &( (rhit)->det()->surface() );
287  LocalVector shitLMom;
288  LocalPoint shitLPos;
289  if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())) {
290  double rechitmatchedx = rhit->localPosition().x();
291  double rechitmatchedy = rhit->localPosition().y();
292  double mindist = 999999;
293  float distx, disty;
294  std::pair<LocalPoint,LocalVector> closestPair;
295  const StripGeomDetUnit* stripDet =(StripGeomDetUnit*) ((const GluedGeomDet *)(rhit)->det())->stereoDet();
296  const BoundPlane& plane = (rhit)->det()->surface();
297  for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++) {
298  //project simhit;
299  std::pair<LocalPoint,LocalVector> hitPair = projectHit((*m),stripDet,plane);
300  distx = fabs(rechitmatchedx - hitPair.first.x());
301  disty = fabs(rechitmatchedy - hitPair.first.y());
302  double dist = distx*distx+disty*disty;
303  if(sqrt(dist)<mindist){
304  mindist = dist;
305  closestPair = hitPair;
306  }
307  }
308  shitLPos = closestPair.first;
309  shitLMom = closestPair.second;
310  } else {
311  shitLPos = shit.localPosition();
312  shitLMom = shit.momentumAtEntry();
313  }
314  GlobalVector shitGMom = surf->toGlobal(shitLMom);
315  GlobalPoint shitGPos = surf->toGlobal(shitLPos);
316 // if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
317 // double rechitmatchedx = rhit->localPosition().x();
318 // double rechitmatchedy = rhit->localPosition().y();
319 // double mindist = 999999;
320 // double distx, disty;
321 // std::pair<LocalPoint,LocalVector> closestPair;
322 // const StripGeomDetUnit* stripDet =(StripGeomDetUnit*) ((const GluedGeomDet *)(rhit)->det())->stereoDet();
323 // const BoundPlane& plane = (rhit)->det()->surface();
324 // for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
325 // const PSimHit& hit = *m;
326 // const StripTopology& topol = stripDet->specificTopology();
327 // GlobalPoint globalpos = stripDet->surface().toGlobal(hit.localPosition());
328 // LocalPoint localHit = plane.toLocal(globalpos);
329 // //track direction
330 // LocalVector locdir=hit.localDirection();
331 // //rotate track in new frame
332 // GlobalVector globaldir= stripDet->surface().toGlobal(locdir);
333 // LocalVector dir=plane.toLocal(globaldir);
334 // float scale = -localHit.z() / dir.z();
335 // LocalPoint projectedPos = localHit + scale*dir;
336 // float selfAngle = topol.stripAngle( topol.strip( hit.localPosition()));
337 // LocalVector stripDir( sin(selfAngle), cos(selfAngle), 0); // vector along strip in hit frame
338 // LocalVector localStripDir = LocalVector( plane.toLocal(stripDet->surface().toGlobal( stripDir)));
339 // std::pair<LocalPoint,LocalVector> hitPair( projectedPos, localStripDir);
340 // distx = fabs(rechitmatchedx - hitPair.first.x());
341 // disty = fabs(rechitmatchedy - hitPair.first.y());
342 // double dist = distx*distx+disty*disty;
343 // if(sqrt(dist)<mindist){
344 // mindist = dist;
345 // closestPair = hitPair;
346 // }
347 // }
348 // shitLPos = closestPair.first;
349 // shitLMom = closestPair.second;
350 // } else {
351 // shitLPos = shit.localPosition();
352 // shitLMom = shit.momentumAtEntry();
353 // }
354 // GlobalVector shitGMom = surf->toGlobal(shitLMom);
355 // GlobalPoint shitGPos = surf->toGlobal(shitLPos);
356 
357  GlobalVector tsosGMom = currentState.globalMomentum();
358  GlobalError tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
359  GlobalPoint tsosGPos = currentState.globalPosition();
360  GlobalError tsosGPEr = currentState.cartesianError().position();
361 
362  GlobalPoint rhitGPos = (rhit)->globalPosition();
363  GlobalError rhitGPEr = (rhit)->globalPositionError();
364 
365  double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/sqrt(rhitGPEr.cxx());
366  double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/sqrt(rhitGPEr.cyy());
367  double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/sqrt(rhitGPEr.czz());
368  //double pullGPX_rs = (rhitGPos.x()-shitGPos.x());
369  //double pullGPY_rs = (rhitGPos.y()-shitGPos.y());
370  //double pullGPZ_rs = (rhitGPos.z()-shitGPos.z());
371 
372  LogTrace("TestSmoothHits") << "rs" << std::endl;
373 
374  title.str("");
375  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
376  hPullGP_X_rs[title.str()]->Fill( pullGPX_rs );
377  title.str("");
378  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
379  hPullGP_Y_rs[title.str()]->Fill( pullGPY_rs );
380  title.str("");
381  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
382  hPullGP_Z_rs[title.str()]->Fill( pullGPZ_rs );
383 
384  double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx());
385  double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy());
386  double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/sqrt(tsosGPEr.czz()+rhitGPEr.czz());
387  //double pullGPX_tr = (tsosGPos.x()-rhitGPos.x());
388  //double pullGPY_tr = (tsosGPos.y()-rhitGPos.y());
389  //double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z());
390 
391  LogTrace("TestSmoothHits") << "tr" << std::endl;
392 
393  title.str("");
394  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
395  hPullGP_X_tr[title.str()]->Fill( pullGPX_tr );
396  title.str("");
397  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
398  hPullGP_Y_tr[title.str()]->Fill( pullGPY_tr );
399  title.str("");
400  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
401  hPullGP_Z_tr[title.str()]->Fill( pullGPZ_tr );
402 
403  double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/sqrt(tsosGPEr.cxx());
404  double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/sqrt(tsosGPEr.cyy());
405  double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/sqrt(tsosGPEr.czz());
406  //double pullGPX_ts = (tsosGPos.x()-shitGPos.x());
407  //double pullGPY_ts = (tsosGPos.y()-shitGPos.y());
408  //double pullGPZ_ts = (tsosGPos.z()-shitGPos.z());
409 
410  LogTrace("TestSmoothHits") << "ts1" << std::endl;
411 
412  title.str("");
413  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
414  hPullGP_X_ts[title.str()]->Fill( pullGPX_ts );
415  title.str("");
416  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
417  hPullGP_Y_ts[title.str()]->Fill( pullGPY_ts );
418  title.str("");
419  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
420  hPullGP_Z_ts[title.str()]->Fill( pullGPZ_ts );
421 
422  double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/sqrt(tsosGMEr.cxx());
423  double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/sqrt(tsosGMEr.cyy());
424  double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/sqrt(tsosGMEr.czz());
425  //double pullGMX_ts = (tsosGMom.x()-shitGMom.x());
426  //double pullGMY_ts = (tsosGMom.y()-shitGMom.y());
427  //double pullGMZ_ts = (tsosGMom.z()-shitGMom.z());
428 
429  LogTrace("TestSmoothHits") << "ts2" << std::endl;
430 
431  title.str("");
432  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
433  hPullGM_X_ts[title.str()]->Fill( pullGMX_ts );
434  title.str("");
435  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
436  hPullGM_Y_ts[title.str()]->Fill( pullGMY_ts );
437  title.str("");
438  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
439  hPullGM_Z_ts[title.str()]->Fill( pullGMZ_ts );
440 
441  if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
442  //mono
443  LogTrace("TestSmoothHits") << "MONO HIT" << std::endl;
444  auto m = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->monoHit();
445  CTTRHp tMonoHit =
446  theBuilder->build(&m);
447  if (tMonoHit==0) continue;
448  vector<PSimHit> assMonoSimHits = hitAssociator->associateHit(*tMonoHit->hit());
449  if (assMonoSimHits.size()==0) continue;
450  const PSimHit sMonoHit = *(assSimHits.begin());
451  const Surface * monoSurf = &( tMonoHit->det()->surface() );
452  if (monoSurf==0) continue;
453  TSOS monoState = thePropagator->propagate(lastState,*monoSurf);
454  if (monoState.isValid()==0) continue;
455 
456  LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
457  GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
458  LocalPoint monoShitLPos = sMonoHit.localPosition();
459  GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
460 
461  GlobalVector monoTsosGMom = monoState.globalMomentum();
462  GlobalError monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
463  GlobalPoint monoTsosGPos = monoState.globalPosition();
464  GlobalError monoTsosGPEr = monoState.cartesianError().position();
465 
466  GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
467  GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
468 
469  double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/sqrt(monoRhitGPEr.cxx());
470  double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/sqrt(monoRhitGPEr.cyy());
471  double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/sqrt(monoRhitGPEr.czz());
472  //double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x());
473  //double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y());
474  //double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z());
475 
476  title.str("");
477  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
478  hPullGP_X_rs_mono[title.str()]->Fill( pullGPX_rs_mono );
479  title.str("");
480  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
481  hPullGP_Y_rs_mono[title.str()]->Fill( pullGPY_rs_mono );
482  title.str("");
483  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
484  hPullGP_Z_rs_mono[title.str()]->Fill( pullGPZ_rs_mono );
485 
486  double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx());
487  double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy());
488  double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz());
489  //double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x());
490  //double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y());
491  //double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z());
492 
493  title.str("");
494  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
495  hPullGP_X_tr_mono[title.str()]->Fill( pullGPX_tr_mono );
496  title.str("");
497  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
498  hPullGP_Y_tr_mono[title.str()]->Fill( pullGPY_tr_mono );
499  title.str("");
500  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
501  hPullGP_Z_tr_mono[title.str()]->Fill( pullGPZ_tr_mono );
502 
503  double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/sqrt(monoTsosGPEr.cxx());
504  double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/sqrt(monoTsosGPEr.cyy());
505  double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/sqrt(monoTsosGPEr.czz());
506  //double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x());
507  //double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y());
508  //double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z());
509 
510  title.str("");
511  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
512  hPullGP_X_ts_mono[title.str()]->Fill( pullGPX_ts_mono );
513  title.str("");
514  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
515  hPullGP_Y_ts_mono[title.str()]->Fill( pullGPY_ts_mono );
516  title.str("");
517  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
518  hPullGP_Z_ts_mono[title.str()]->Fill( pullGPZ_ts_mono );
519 
520  double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/sqrt(monoTsosGMEr.cxx());
521  double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/sqrt(monoTsosGMEr.cyy());
522  double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/sqrt(monoTsosGMEr.czz());
523  //double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x());
524  //double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y());
525  //double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z());
526 
527  title.str("");
528  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
529  hPullGM_X_ts_mono[title.str()]->Fill( pullGMX_ts_mono );
530  title.str("");
531  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
532  hPullGM_Y_ts_mono[title.str()]->Fill( pullGMY_ts_mono );
533  title.str("");
534  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
535  hPullGM_Z_ts_mono[title.str()]->Fill( pullGMZ_ts_mono );
536 
537  //stereo
538  LogTrace("TestSmoothHits") << "STEREO HIT" << std::endl;
539  auto s = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->stereoHit();
540  CTTRHp tStereoHit =
541  theBuilder->build(&s);
542  if (tStereoHit==0) continue;
543  vector<PSimHit> assStereoSimHits = hitAssociator->associateHit(*tStereoHit->hit());
544  if (assStereoSimHits.size()==0) continue;
545  const PSimHit sStereoHit = *(assSimHits.begin());
546  const Surface * stereoSurf = &( tStereoHit->det()->surface() );
547  if (stereoSurf==0) continue;
548  TSOS stereoState = thePropagator->propagate(lastState,*stereoSurf);
549  if (stereoState.isValid()==0) continue;
550 
551  LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
552  GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
553  LocalPoint stereoShitLPos = sStereoHit.localPosition();
554  GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
555 
556  GlobalVector stereoTsosGMom = stereoState.globalMomentum();
557  GlobalError stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
558  GlobalPoint stereoTsosGPos = stereoState.globalPosition();
559  GlobalError stereoTsosGPEr = stereoState.cartesianError().position();
560 
561  GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
562  GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
563 
564  double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/sqrt(stereoRhitGPEr.cxx());
565  double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/sqrt(stereoRhitGPEr.cyy());
566  double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/sqrt(stereoRhitGPEr.czz());
567  //double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x());
568  //double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y());
569  //double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z());
570 
571  title.str("");
572  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
573  hPullGP_X_rs_stereo[title.str()]->Fill( pullGPX_rs_stereo );
574  title.str("");
575  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
576  hPullGP_Y_rs_stereo[title.str()]->Fill( pullGPY_rs_stereo );
577  title.str("");
578  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
579  hPullGP_Z_rs_stereo[title.str()]->Fill( pullGPZ_rs_stereo );
580 
581  double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx());
582  double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy());
583  double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz());
584  //double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x());
585  //double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y());
586  //double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z());
587 
588  title.str("");
589  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
590  hPullGP_X_tr_stereo[title.str()]->Fill( pullGPX_tr_stereo );
591  title.str("");
592  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
593  hPullGP_Y_tr_stereo[title.str()]->Fill( pullGPY_tr_stereo );
594  title.str("");
595  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
596  hPullGP_Z_tr_stereo[title.str()]->Fill( pullGPZ_tr_stereo );
597 
598  double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/sqrt(stereoTsosGPEr.cxx());
599  double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/sqrt(stereoTsosGPEr.cyy());
600  double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/sqrt(stereoTsosGPEr.czz());
601  //double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x());
602  //double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y());
603  //double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z());
604 
605  title.str("");
606  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
607  hPullGP_X_ts_stereo[title.str()]->Fill( pullGPX_ts_stereo );
608  title.str("");
609  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
610  hPullGP_Y_ts_stereo[title.str()]->Fill( pullGPY_ts_stereo );
611  title.str("");
612  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
613  hPullGP_Z_ts_stereo[title.str()]->Fill( pullGPZ_ts_stereo );
614 
615  double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/sqrt(stereoTsosGMEr.cxx());
616  double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/sqrt(stereoTsosGMEr.cyy());
617  double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/sqrt(stereoTsosGMEr.czz());
618  //double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x());
619  //double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y());
620  //double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z());
621 
622  title.str("");
623  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
624  hPullGM_X_ts_stereo[title.str()]->Fill( pullGMX_ts_stereo );
625  title.str("");
626  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
627  hPullGM_Y_ts_stereo[title.str()]->Fill( pullGMY_ts_stereo );
628  title.str("");
629  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
630  hPullGM_Z_ts_stereo[title.str()]->Fill( pullGMZ_ts_stereo );
631  }
632  lastState = updatedState;
633  //#endif
634  }
635  }
636  delete hitAssociator;
637  LogTrace("TestSmoothHits") << "end of event" << std::endl;
638 }
639 
641  //file->Write();
642  TDirectory * chi2i = file->mkdir("Chi2_Increment");
643 
644  TDirectory * gp_ts = file->mkdir("GP_TSOS-SimHit");
645  TDirectory * gm_ts = file->mkdir("GM_TSOS-SimHit");
646  TDirectory * gp_tr = file->mkdir("GP_TSOS-RecHit");
647  TDirectory * gp_rs = file->mkdir("GP_RecHit-SimHit");
648 
649  TDirectory * gp_tsx = gp_ts->mkdir("X");
650  TDirectory * gp_tsy = gp_ts->mkdir("Y");
651  TDirectory * gp_tsz = gp_ts->mkdir("Z");
652  TDirectory * gm_tsx = gm_ts->mkdir("X");
653  TDirectory * gm_tsy = gm_ts->mkdir("Y");
654  TDirectory * gm_tsz = gm_ts->mkdir("Z");
655  TDirectory * gp_trx = gp_tr->mkdir("X");
656  TDirectory * gp_try = gp_tr->mkdir("Y");
657  TDirectory * gp_trz = gp_tr->mkdir("Z");
658  TDirectory * gp_rsx = gp_rs->mkdir("X");
659  TDirectory * gp_rsy = gp_rs->mkdir("Y");
660  TDirectory * gp_rsz = gp_rs->mkdir("Z");
661 
662  TDirectory * gp_tsx_mono = gp_ts->mkdir("MONOX");
663  TDirectory * gp_tsy_mono = gp_ts->mkdir("MONOY");
664  TDirectory * gp_tsz_mono = gp_ts->mkdir("MONOZ");
665  TDirectory * gm_tsx_mono = gm_ts->mkdir("MONOX");
666  TDirectory * gm_tsy_mono = gm_ts->mkdir("MONOY");
667  TDirectory * gm_tsz_mono = gm_ts->mkdir("MONOZ");
668  TDirectory * gp_trx_mono = gp_tr->mkdir("MONOX");
669  TDirectory * gp_try_mono = gp_tr->mkdir("MONOY");
670  TDirectory * gp_trz_mono = gp_tr->mkdir("MONOZ");
671  TDirectory * gp_rsx_mono = gp_rs->mkdir("MONOX");
672  TDirectory * gp_rsy_mono = gp_rs->mkdir("MONOY");
673  TDirectory * gp_rsz_mono = gp_rs->mkdir("MONOZ");
674 
675  TDirectory * gp_tsx_stereo = gp_ts->mkdir("STEREOX");
676  TDirectory * gp_tsy_stereo = gp_ts->mkdir("STEREOY");
677  TDirectory * gp_tsz_stereo = gp_ts->mkdir("STEREOZ");
678  TDirectory * gm_tsx_stereo = gm_ts->mkdir("STEREOX");
679  TDirectory * gm_tsy_stereo = gm_ts->mkdir("STEREOY");
680  TDirectory * gm_tsz_stereo = gm_ts->mkdir("STEREOZ");
681  TDirectory * gp_trx_stereo = gp_tr->mkdir("STEREOX");
682  TDirectory * gp_try_stereo = gp_tr->mkdir("STEREOY");
683  TDirectory * gp_trz_stereo = gp_tr->mkdir("STEREOZ");
684  TDirectory * gp_rsx_stereo = gp_rs->mkdir("STEREOX");
685  TDirectory * gp_rsy_stereo = gp_rs->mkdir("STEREOY");
686  TDirectory * gp_rsz_stereo = gp_rs->mkdir("STEREOZ");
687 
688  chi2i->cd();
689  hTotChi2Increment->Write();
690  hChi2_vs_Process->Write();
691  hChi2_vs_clsize->Write();
692  for (int i=0; i!=6; i++)
693  for (int j=0; j!=9; j++){
694  if (i==0 && j>2) break;
695  if (i==1 && j>1) break;
696  if (i==2 && j>3) break;
697  if (i==3 && j>2) break;
698  if (i==4 && j>5) break;
699  if (i==5 && j>8) break;
700  chi2i->cd();
701  title.str("");
702  title << "Chi2Increment_" << i+1 << "-" << j+1;
703  hChi2Increment[title.str()]->Write();
704 
705  gp_ts->cd();
706  gp_tsx->cd();
707  title.str("");
708  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
709  hPullGP_X_ts[title.str()]->Write();
710  gp_tsy->cd();
711  title.str("");
712  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts";
713  hPullGP_Y_ts[title.str()]->Write();
714  gp_tsz->cd();
715  title.str("");
716  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts";
717  hPullGP_Z_ts[title.str()]->Write();
718 
719  gm_ts->cd();
720  gm_tsx->cd();
721  title.str("");
722  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
723  hPullGM_X_ts[title.str()]->Write();
724  gm_tsy->cd();
725  title.str("");
726  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
727  hPullGM_Y_ts[title.str()]->Write();
728  gm_tsz->cd();
729  title.str("");
730  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
731  hPullGM_Z_ts[title.str()]->Write();
732 
733  gp_tr->cd();
734  gp_trx->cd();
735  title.str("");
736  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
737  hPullGP_X_tr[title.str()]->Write();
738  gp_try->cd();
739  title.str("");
740  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
741  hPullGP_Y_tr[title.str()]->Write();
742  gp_trz->cd();
743  title.str("");
744  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
745  hPullGP_Z_tr[title.str()]->Write();
746 
747  gp_rs->cd();
748  gp_rsx->cd();
749  title.str("");
750  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
751  hPullGP_X_rs[title.str()]->Write();
752  gp_rsy->cd();
753  title.str("");
754  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
755  hPullGP_Y_rs[title.str()]->Write();
756  gp_rsz->cd();
757  title.str("");
758  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
759  hPullGP_Z_rs[title.str()]->Write();
760 
761  if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
762  //mono
763  gp_ts->cd();
764  gp_tsx_mono->cd();
765  title.str("");
766  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
767  hPullGP_X_ts_mono[title.str()]->Write();
768  gp_tsy_mono->cd();
769  title.str("");
770  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
771  hPullGP_Y_ts_mono[title.str()]->Write();
772  gp_tsz_mono->cd();
773  title.str("");
774  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
775  hPullGP_Z_ts_mono[title.str()]->Write();
776 
777  gm_ts->cd();
778  gm_tsx_mono->cd();
779  title.str("");
780  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
781  hPullGM_X_ts_mono[title.str()]->Write();
782  gm_tsy_mono->cd();
783  title.str("");
784  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
785  hPullGM_Y_ts_mono[title.str()]->Write();
786  gm_tsz_mono->cd();
787  title.str("");
788  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
789  hPullGM_Z_ts_mono[title.str()]->Write();
790 
791  gp_tr->cd();
792  gp_trx_mono->cd();
793  title.str("");
794  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
795  hPullGP_X_tr_mono[title.str()]->Write();
796  gp_try_mono->cd();
797  title.str("");
798  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
799  hPullGP_Y_tr_mono[title.str()]->Write();
800  gp_trz_mono->cd();
801  title.str("");
802  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
803  hPullGP_Z_tr_mono[title.str()]->Write();
804 
805  gp_rs->cd();
806  gp_rsx_mono->cd();
807  title.str("");
808  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
809  hPullGP_X_rs_mono[title.str()]->Write();
810  gp_rsy_mono->cd();
811  title.str("");
812  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
813  hPullGP_Y_rs_mono[title.str()]->Write();
814  gp_rsz_mono->cd();
815  title.str("");
816  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
817  hPullGP_Z_rs_mono[title.str()]->Write();
818 
819  //stereo
820  gp_ts->cd();
821  gp_tsx_stereo->cd();
822  title.str("");
823  title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
824  hPullGP_X_ts_stereo[title.str()]->Write();
825  gp_tsy_stereo->cd();
826  title.str("");
827  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
828  hPullGP_Y_ts_stereo[title.str()]->Write();
829  gp_tsz_stereo->cd();
830  title.str("");
831  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
832  hPullGP_Z_ts_stereo[title.str()]->Write();
833 
834  gm_ts->cd();
835  gm_tsx_stereo->cd();
836  title.str("");
837  title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
838  hPullGM_X_ts_stereo[title.str()]->Write();
839  gm_tsy_stereo->cd();
840  title.str("");
841  title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
842  hPullGM_Y_ts_stereo[title.str()]->Write();
843  gm_tsz_stereo->cd();
844  title.str("");
845  title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
846  hPullGM_Z_ts_stereo[title.str()]->Write();
847 
848  gp_tr->cd();
849  gp_trx_stereo->cd();
850  title.str("");
851  title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
852  hPullGP_X_tr_stereo[title.str()]->Write();
853  gp_try_stereo->cd();
854  title.str("");
855  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
856  hPullGP_Y_tr_stereo[title.str()]->Write();
857  gp_trz_stereo->cd();
858  title.str("");
859  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
860  hPullGP_Z_tr_stereo[title.str()]->Write();
861 
862  gp_rs->cd();
863  gp_rsx_stereo->cd();
864  title.str("");
865  title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
866  hPullGP_X_rs_stereo[title.str()]->Write();
867  gp_rsy_stereo->cd();
868  title.str("");
869  title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
870  hPullGP_Y_rs_stereo[title.str()]->Write();
871  gp_rsz_stereo->cd();
872  title.str("");
873  title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
874  hPullGP_Z_rs_stereo[title.str()]->Write();
875  }
876  }
877 
878  file->Close();
879 }
880 
881 //needed by to do the residual for matched hits
882 //taken from SiStripTrackingRecHitsValid.cc
883 std::pair<LocalPoint,LocalVector>
884 TestSmoothHits::projectHit( const PSimHit& hit, const StripGeomDetUnit* stripDet, const BoundPlane& plane)
885 {
886  const StripTopology& topol = stripDet->specificTopology();
887  GlobalPoint globalpos= stripDet->surface().toGlobal(hit.localPosition());
888  LocalPoint localHit = plane.toLocal(globalpos);
889  //track direction
890  LocalVector locdir=hit.localDirection();
891  //rotate track in new frame
892 
893  GlobalVector globaldir= stripDet->surface().toGlobal(locdir);
894  LocalVector dir=plane.toLocal(globaldir);
895  float scale = -localHit.z() / dir.z();
896 
897  LocalPoint projectedPos = localHit + scale*dir;
898 
899  float selfAngle = topol.stripAngle( topol.strip( hit.localPosition()));
900 
901  LocalVector stripDir( sin(selfAngle), cos(selfAngle), 0); // vector along strip in hit frame
902 
903  LocalVector localStripDir( plane.toLocal(stripDet->surface().toGlobal( stripDir)));
904 
905  return std::pair<LocalPoint,LocalVector>( projectedPos, localStripDir);
906 }
907 
910 
std::map< std::string, TH1F * > hPullGP_Z_rs_stereo
dbl * delta
Definition: mlp_gen.cc:36
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:78
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
std::string srcName
TrackerHitAssociator * hitAssociator
std::map< std::string, TH1F * > hPullGP_X_ts_mono
std::map< std::string, TH1F * > hPullGM_Y_ts_mono
edm::Handle< TrackCandidateCollection > theTCCollection
virtual float stripAngle(float strip) const =0
std::stringstream title
range recHits() const
TestSmoothHits(const edm::ParameterSet &)
LocalVector momentumAtEntry() const
The momentum of the track that produced the hit, at entry point.
Definition: PSimHit.h:47
std::map< std::string, TH1F * > hChi2Increment
virtual void endJob()
std::map< std::string, TH1F * > hPullGP_X_rs_mono
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
std::string propagatorName
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
TrajectorySeed const & seed() const
std::map< std::string, TH1F * > hPullGM_Z_ts_stereo
const CartesianTrajectoryError cartesianError() const
T y() const
Definition: PV3DBase.h:62
GlobalPoint globalPosition() const
std::map< std::string, TH1F * > hPullGP_Y_rs
std::map< std::string, TH1F * > hPullGM_X_ts_mono
std::map< std::string, TH1F * > hPullGP_Z_tr
std::map< std::string, TH1F * > hPullGP_X_ts
std::pair< const_iterator, const_iterator > range
std::string fname
std::map< std::string, TH1F * > hPullGP_Y_tr_stereo
std::map< std::string, TH1F * > hPullGP_Z_ts
virtual const StripTopology & specificTopology() const
Returns a reference to the strip proxy topology.
std::map< std::string, TH1F * > hPullGP_X_rs_stereo
std::vector< ConstRecHitPointer > RecHitContainer
edm::ESHandle< TrajectorySmoother > smooth
std::map< std::string, TH1F * > hPullGP_Y_rs_mono
virtual float strip(const LocalPoint &) const =0
std::map< std::string, TH1F * > hPullGM_Y_ts_stereo
edm::ESHandle< TrajectoryFitter > fit
PTrajectoryStateOnDet const & trajectoryStateOnDet() const
TH2F * hChi2_vs_Process
std::map< std::string, TH1F * > hPullGP_Z_rs_mono
edm::ESHandle< MagneticField > theMF
TH2F * hChi2_vs_clsize
int iEvent
Definition: GenABIO.cc:243
std::map< std::string, TH1F * > hPullGP_Y_ts_mono
std::map< std::string, TH1F * > hPullGM_Y_ts
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
std::map< std::string, TH1F * > hPullGP_Z_tr_mono
Local3DPoint localPosition() const
Definition: PSimHit.h:44
T sqrt(T t)
Definition: SSEVec.h:46
LocalPoint toLocal(const GlobalPoint &gp) const
T z() const
Definition: PV3DBase.h:63
tuple result
Definition: query.py:137
std::map< std::string, TH1F * > hPullGM_Z_ts
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
std::map< std::string, TH1F * > hPullGP_X_tr
int j
Definition: DBlmapReader.cc:9
const edm::ParameterSet conf_
TransientTrackingRecHit::ConstRecHitPointer CTTRHp
Definition: CkfDebugger.h:44
unsigned int detId() const
std::map< std::string, TH1F * > hPullGP_Y_rs_stereo
std::string sname
virtual void analyze(const edm::Event &, const edm::EventSetup &)
edm::ESHandle< TransientTrackingRecHitBuilder > theBuilder
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
virtual void beginRun(edm::Run &run, const edm::EventSetup &)
#define LogTrace(id)
std::map< std::string, TH1F * > hPullGP_Y_tr_mono
const AlgebraicSymMatrix66 & matrix() const
std::map< std::string, TH1F * > hPullGP_X_ts_stereo
std::map< std::string, TH1F * > hPullGP_Y_ts
LocalVector localDirection() const
Obsolete. Same as momentumAtEntry().unit(), for backward compatibility.
Definition: PSimHit.h:52
std::map< std::string, TH1F * > hPullGP_Z_ts_stereo
Definition: DetId.h:20
TrajectoryStateOnSurface TSOS
Definition: TestHits.cc:23
std::map< std::string, TH1F * > hPullGM_Z_ts_mono
TrajectoryStateOnSurface transientState(const PTrajectoryStateOnDet &ts, const Surface *surface, const MagneticField *field)
const GlobalError position() const
Position error submatrix.
std::map< std::string, TH1F * > hPullGP_X_rs
TH1F * hTotChi2Increment
edm::ESHandle< Propagator > thePropagator
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:62
std::map< std::string, TH1F * > hPullGP_X_tr_stereo
unsigned short processType() const
Definition: PSimHit.h:118
std::string builderName
char state
Definition: procUtils.cc:75
edm::ESHandle< TrackerGeometry > theG
std::map< std::string, TH1F * > hPullGP_X_tr_mono
std::pair< LocalPoint, LocalVector > projectHit(const PSimHit &, const StripGeomDetUnit *, const BoundPlane &)
GlobalVector globalMomentum() const
const BoundPlane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:35
std::vector< PSimHit > associateHit(const TrackingRecHit &thit)
std::map< std::string, TH1F * > hPullGP_Z_tr_stereo
std::map< std::string, TH1F * > hPullGP_Z_rs
dbl *** dir
Definition: mlp_gen.cc:35
std::map< std::string, TH1F * > hPullGP_Y_ts_stereo
std::map< std::string, TH1F * > hPullGP_Z_ts_mono
std::map< std::string, TH1F * > hPullGM_X_ts
T x() const
Definition: PV3DBase.h:61
std::map< std::string, TH1F * > hPullGM_X_ts_stereo
std::map< std::string, TH1F * > hPullGP_Y_tr
Definition: Run.h:33
Pixel Reconstructed Hit.