CMS 3D CMS Logo

TestHits.cc
Go to the documentation of this file.
1 #include "TestHits.h"
2 
11 #include <TDirectory.h>
14 
18 
21 using namespace std;
22 using namespace edm;
23 
24 TestHits::TestHits(const edm::ParameterSet& iConfig) : trackerHitAssociatorConfig_(consumesCollector()) {
25  LogTrace("TestHits") << iConfig << std::endl;
26  propagatorName = iConfig.getParameter<std::string>("Propagator");
27  builderName = iConfig.getParameter<std::string>("TTRHBuilder");
28  srcName = iConfig.getParameter<std::string>("src");
29  fname = iConfig.getParameter<std::string>("Fitter");
30  mineta = iConfig.getParameter<double>("mineta");
31  maxeta = iConfig.getParameter<double>("maxeta");
32 
33  theGToken = esConsumes<edm::Transition::BeginRun>();
34  theMFToken = esConsumes<edm::Transition::BeginRun>();
35  thePropagatorToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", propagatorName));
36  theBuilderToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", builderName));
37  fitToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", fname));
39 
41 }
42 
44 
45 void TestHits::beginRun(edm::Run const& run, const edm::EventSetup& iSetup) {
46  theG = iSetup.getHandle(theGToken);
47  theMF = iSetup.getHandle(theMFToken);
50  fit = iSetup.getHandle(fitToken);
51 
52  file = new TFile("testhits.root", "recreate");
53  for (int i = 0; i != 6; i++)
54  for (int j = 0; j != 9; j++) {
55  if (i == 0 && j > 2)
56  break;
57  if (i == 1 && j > 1)
58  break;
59  if (i == 2 && j > 3)
60  break;
61  if (i == 3 && j > 2)
62  break;
63  if (i == 4 && j > 5)
64  break;
65  if (i == 5 && j > 8)
66  break;
67  title.str("");
68  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts";
69  hPullGP_X_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
70  title.str("");
71  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts";
72  hPullGP_Y_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
73  title.str("");
74  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts";
75  hPullGP_Z_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
76  title.str("");
77  title << "Chi2Increment_" << i + 1 << "-" << j + 1;
78  hChi2Increment[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, 0, 100);
79 
80  title.str("");
81  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts";
82  hPullGM_X_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
83  title.str("");
84  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts";
85  hPullGM_Y_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
86  title.str("");
87  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts";
88  hPullGM_Z_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
89 
90  title.str("");
91  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr";
92  hPullGP_X_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
93  title.str("");
94  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr";
95  hPullGP_Y_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
96  title.str("");
97  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr";
98  hPullGP_Z_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
99 
100  title.str("");
101  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs";
102  hPullGP_X_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
103  title.str("");
104  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs";
105  hPullGP_Y_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
106  title.str("");
107  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs";
108  hPullGP_Z_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
109 
110  if (((i == 2 || i == 4) && (j == 0 || j == 1)) || (i == 3 || i == 5)) {
111  //mono
112  title.str("");
113  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
114  hPullGP_X_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
115  title.str("");
116  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
117  hPullGP_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
118  title.str("");
119  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
120  hPullGP_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
121 
122  title.str("");
123  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
124  hPullGM_X_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
125  title.str("");
126  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
127  hPullGM_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
128  title.str("");
129  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
130  hPullGM_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
131 
132  title.str("");
133  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_mono";
134  hPullGP_X_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
135  title.str("");
136  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_mono";
137  hPullGP_Y_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
138  title.str("");
139  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_mono";
140  hPullGP_Z_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
141 
142  title.str("");
143  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_mono";
144  hPullGP_X_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
145  title.str("");
146  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_mono";
147  hPullGP_Y_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
148  title.str("");
149  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_mono";
150  hPullGP_Z_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
151 
152  //stereo
153  title.str("");
154  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
155  hPullGP_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
156  title.str("");
157  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
158  hPullGP_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
159  title.str("");
160  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
161  hPullGP_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
162 
163  title.str("");
164  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
165  hPullGM_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
166  title.str("");
167  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
168  hPullGM_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
169  title.str("");
170  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
171  hPullGM_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
172 
173  title.str("");
174  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_stereo";
175  hPullGP_X_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
176  title.str("");
177  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_stereo";
178  hPullGP_Y_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
179  title.str("");
180  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_stereo";
181  hPullGP_Z_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
182 
183  title.str("");
184  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_stereo";
185  hPullGP_X_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
186  title.str("");
187  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_stereo";
188  hPullGP_Y_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
189  title.str("");
190  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_stereo";
191  hPullGP_Z_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
192  }
193  }
194  hTotChi2Increment = new TH1F("TotChi2Increment", "TotChi2Increment", 1000, 0, 100);
195  hProcess_vs_Chi2 = new TH2F("Process_vs_Chi2", "Process_vs_Chi2", 1000, 0, 100, 17, -0.5, 16.5);
196  hClsize_vs_Chi2 = new TH2F("Clsize_vs_Chi2", "Clsize_vs_Chi2", 1000, 0, 100, 17, -0.5, 16.5);
197 }
198 
199 void TestHits::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
200  //Retrieve tracker topology from geometry
202 
203  LogTrace("TestHits") << "\nnew event";
204 
207 
209 
210  for (TrackCandidateCollection::const_iterator i = theTCCollection->begin(); i != theTCCollection->end(); i++) {
211  LogTrace("TestHits") << "\n*****************new candidate*****************" << std::endl;
212 
213  const TrackCandidate* theTC = &(*i);
215 
216  //convert PTrajectoryStateOnDet to TrajectoryStateOnSurface
217 
218  DetId detId(state.detId());
219  TrajectoryStateOnSurface theTSOS =
221 
222  if (theTSOS.globalMomentum().eta() > maxeta || theTSOS.globalMomentum().eta() < mineta)
223  continue;
224 
225  //convert the TrackingRecHit vector to a TransientTrackingRecHit vector
227 
228  for (auto const& recHit : theTC->recHits()) {
229  hits.push_back(theBuilder->build(&recHit));
230  }
231 
232  //call the fitter
233  std::vector<Trajectory> result = fit->fit(theTC->seed(), hits, theTSOS);
234  if (result.empty())
235  continue;
236  std::vector<TrajectoryMeasurement> vtm = result[0].measurements();
237  double tchi2 = 0;
238 
239  TSOS lastState = theTSOS;
240  for (std::vector<TrajectoryMeasurement>::iterator tm = vtm.begin(); tm != vtm.end(); tm++) {
242  if ((rhit)->isValid() == 0 && rhit->det() != nullptr)
243  continue;
244  LogTrace("TestHits") << "*****************new hit*****************";
245 
246  int subdetId = rhit->det()->geographicalId().subdetId();
247  DetId id = rhit->det()->geographicalId();
248  int layerId = tTopo->layer(id);
249  LogTrace("TestHits") << "subdetId=" << subdetId << " layerId=" << layerId;
250 
251  double delta = 99999;
252  LocalPoint rhitLPv = rhit->localPosition();
253 
254  std::vector<PSimHit> assSimHits = hitAssociator.associateHit(*(rhit->hit()));
255  if (assSimHits.empty())
256  continue;
257  PSimHit shit;
258  for (std::vector<PSimHit>::const_iterator m = assSimHits.begin(); m < assSimHits.end(); m++) {
259  if ((m->localPosition() - rhitLPv).mag() < delta) {
260  shit = *m;
261  delta = (m->localPosition() - rhitLPv).mag();
262  }
263  }
264 
265  TSOS currentState = tm->forwardPredictedState();
266  if (tm->backwardPredictedState().isValid())
267  currentState = combiner(tm->backwardPredictedState(), tm->forwardPredictedState());
268  TSOS updatedState = tm->updatedState();
269  tchi2 += tm->estimate();
270 
271  //plot chi2 increment
272  double chi2increment = tm->estimate();
273  LogTrace("TestHits") << "tm->estimate()=" << tm->estimate();
274  title.str("");
275  title << "Chi2Increment_" << subdetId << "-" << layerId;
276  hChi2Increment[title.str()]->Fill(chi2increment);
277  hTotChi2Increment->Fill(chi2increment);
278  hProcess_vs_Chi2->Fill(chi2increment, shit.processType());
279  if (dynamic_cast<const SiPixelRecHit*>(rhit->hit()))
280  hClsize_vs_Chi2->Fill(chi2increment, ((const SiPixelRecHit*)(rhit->hit()))->cluster()->size());
281  if (dynamic_cast<const SiStripRecHit2D*>(rhit->hit()))
282  hClsize_vs_Chi2->Fill(chi2increment, ((const SiStripRecHit2D*)(rhit->hit()))->cluster()->amplitudes().size());
283 
284  //test hits
285  const Surface* surf = &((rhit)->det()->surface());
286  LocalVector shitLMom;
287  LocalPoint shitLPos;
288  if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
289  double rechitmatchedx = rhit->localPosition().x();
290  double rechitmatchedy = rhit->localPosition().y();
291  double mindist = 999999;
292  double distx, disty;
293  std::pair<LocalPoint, LocalVector> closestPair;
294  const StripGeomDetUnit* stripDet = (StripGeomDetUnit*)((const GluedGeomDet*)(rhit)->det())->stereoDet();
295  const BoundPlane& plane = (rhit)->det()->surface();
296  for (std::vector<PSimHit>::const_iterator m = assSimHits.begin(); m < assSimHits.end(); m++) {
297  //project simhit;
298  std::pair<LocalPoint, LocalVector> hitPair = projectHit((*m), stripDet, plane);
299  distx = fabs(rechitmatchedx - hitPair.first.x());
300  disty = fabs(rechitmatchedy - hitPair.first.y());
301  double dist = distx * distx + disty * disty;
302  if (sqrt(dist) < mindist) {
303  mindist = dist;
304  closestPair = hitPair;
305  }
306  }
307  shitLPos = closestPair.first;
308  shitLMom = closestPair.second;
309  } else {
310  shitLPos = shit.localPosition();
311  shitLMom = shit.momentumAtEntry();
312  }
313  GlobalVector shitGMom = surf->toGlobal(shitLMom);
314  GlobalPoint shitGPos = surf->toGlobal(shitLPos);
315 
316  GlobalVector tsosGMom = currentState.globalMomentum();
317  GlobalError tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
318  GlobalPoint tsosGPos = currentState.globalPosition();
319  GlobalError tsosGPEr = currentState.cartesianError().position();
320 
321  GlobalPoint rhitGPos = (rhit)->globalPosition();
322  GlobalError rhitGPEr = (rhit)->globalPositionError();
323 
324  double pullGPX_rs = (rhitGPos.x() - shitGPos.x()) / sqrt(rhitGPEr.cxx());
325  double pullGPY_rs = (rhitGPos.y() - shitGPos.y()) / sqrt(rhitGPEr.cyy());
326  double pullGPZ_rs = (rhitGPos.z() - shitGPos.z()) / sqrt(rhitGPEr.czz());
327  //double pullGPX_rs = (rhitGPos.x()-shitGPos.x());
328  //double pullGPY_rs = (rhitGPos.y()-shitGPos.y());
329  //double pullGPZ_rs = (rhitGPos.z()-shitGPos.z());
330 
331  LogTrace("TestHits") << "rs" << std::endl;
332  LogVerbatim("TestHits") << "assSimHits.size()=" << assSimHits.size();
333  LogVerbatim("TestHits") << "tsos globalPos =" << tsosGPos;
334  LogVerbatim("TestHits") << "sim hit globalPos=" << shitGPos;
335  LogVerbatim("TestHits") << "rec hit globalPos=" << rhitGPos;
336  LogVerbatim("TestHits") << "geographicalId =" << rhit->det()->geographicalId().rawId();
337  LogVerbatim("TestHits") << "surface position =" << surf->position();
338 
339  title.str("");
340  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
341  hPullGP_X_rs[title.str()]->Fill(pullGPX_rs);
342  title.str("");
343  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
344  hPullGP_Y_rs[title.str()]->Fill(pullGPY_rs);
345  title.str("");
346  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
347  hPullGP_Z_rs[title.str()]->Fill(pullGPZ_rs);
348 
349  double pullGPX_tr = (tsosGPos.x() - rhitGPos.x()) / sqrt(tsosGPEr.cxx() + rhitGPEr.cxx());
350  double pullGPY_tr = (tsosGPos.y() - rhitGPos.y()) / sqrt(tsosGPEr.cyy() + rhitGPEr.cyy());
351  double pullGPZ_tr = (tsosGPos.z() - rhitGPos.z()) / sqrt(tsosGPEr.czz() + rhitGPEr.czz());
352  //double pullGPX_tr = (tsosGPos.x()-rhitGPos.x());
353  //double pullGPY_tr = (tsosGPos.y()-rhitGPos.y());
354  //double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z());
355 
356  LogTrace("TestHits") << "tr" << std::endl;
357 
358  title.str("");
359  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
360  hPullGP_X_tr[title.str()]->Fill(pullGPX_tr);
361  title.str("");
362  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
363  hPullGP_Y_tr[title.str()]->Fill(pullGPY_tr);
364  title.str("");
365  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
366  hPullGP_Z_tr[title.str()]->Fill(pullGPZ_tr);
367 
368  double pullGPX_ts = (tsosGPos.x() - shitGPos.x()) / sqrt(tsosGPEr.cxx());
369  double pullGPY_ts = (tsosGPos.y() - shitGPos.y()) / sqrt(tsosGPEr.cyy());
370  double pullGPZ_ts = (tsosGPos.z() - shitGPos.z()) / sqrt(tsosGPEr.czz());
371  //double pullGPX_ts = (tsosGPos.x()-shitGPos.x());
372  //double pullGPY_ts = (tsosGPos.y()-shitGPos.y());
373  //double pullGPZ_ts = (tsosGPos.z()-shitGPos.z());
374 
375  LogTrace("TestHits") << "ts1" << std::endl;
376 
377  title.str("");
378  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
379  hPullGP_X_ts[title.str()]->Fill(pullGPX_ts);
380  title.str("");
381  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
382  hPullGP_Y_ts[title.str()]->Fill(pullGPY_ts);
383  title.str("");
384  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
385  hPullGP_Z_ts[title.str()]->Fill(pullGPZ_ts);
386 
387  double pullGMX_ts = (tsosGMom.x() - shitGMom.x()) / sqrt(tsosGMEr.cxx());
388  double pullGMY_ts = (tsosGMom.y() - shitGMom.y()) / sqrt(tsosGMEr.cyy());
389  double pullGMZ_ts = (tsosGMom.z() - shitGMom.z()) / sqrt(tsosGMEr.czz());
390  //double pullGMX_ts = (tsosGMom.x()-shitGMom.x());
391  //double pullGMY_ts = (tsosGMom.y()-shitGMom.y());
392  //double pullGMZ_ts = (tsosGMom.z()-shitGMom.z());
393 
394  LogTrace("TestHits") << "ts2" << std::endl;
395 
396  title.str("");
397  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
398  hPullGM_X_ts[title.str()]->Fill(pullGMX_ts);
399  title.str("");
400  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
401  hPullGM_Y_ts[title.str()]->Fill(pullGMY_ts);
402  title.str("");
403  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
404  hPullGM_Z_ts[title.str()]->Fill(pullGMZ_ts);
405 
406  if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
407  //mono
408  LogTrace("TestHits") << "MONO HIT" << std::endl;
409  auto m = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->monoHit();
410  CTTRHp tMonoHit = theBuilder->build(&m);
411  if (tMonoHit == nullptr)
412  continue;
413  vector<PSimHit> assMonoSimHits = hitAssociator.associateHit(*tMonoHit->hit());
414  if (assMonoSimHits.empty())
415  continue;
416  const PSimHit sMonoHit = *(assSimHits.begin());
417  const Surface* monoSurf = &(tMonoHit->det()->surface());
418  if (monoSurf == nullptr)
419  continue;
420  TSOS monoState = thePropagator->propagate(lastState, *monoSurf);
421  if (monoState.isValid() == 0)
422  continue;
423 
424  LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
425  GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
426  LocalPoint monoShitLPos = sMonoHit.localPosition();
427  GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
428 
429  GlobalVector monoTsosGMom = monoState.globalMomentum();
430  GlobalError monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
431  GlobalPoint monoTsosGPos = monoState.globalPosition();
432  GlobalError monoTsosGPEr = monoState.cartesianError().position();
433 
434  GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
435  GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
436 
437  double pullGPX_rs_mono = (monoRhitGPos.x() - monoShitGPos.x()) / sqrt(monoRhitGPEr.cxx());
438  double pullGPY_rs_mono = (monoRhitGPos.y() - monoShitGPos.y()) / sqrt(monoRhitGPEr.cyy());
439  double pullGPZ_rs_mono = (monoRhitGPos.z() - monoShitGPos.z()) / sqrt(monoRhitGPEr.czz());
440  //double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x());
441  //double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y());
442  //double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z());
443 
444  title.str("");
445  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
446  hPullGP_X_rs_mono[title.str()]->Fill(pullGPX_rs_mono);
447  title.str("");
448  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
449  hPullGP_Y_rs_mono[title.str()]->Fill(pullGPY_rs_mono);
450  title.str("");
451  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
452  hPullGP_Z_rs_mono[title.str()]->Fill(pullGPZ_rs_mono);
453 
454  double pullGPX_tr_mono = (monoTsosGPos.x() - monoRhitGPos.x()) / sqrt(monoTsosGPEr.cxx() + monoRhitGPEr.cxx());
455  double pullGPY_tr_mono = (monoTsosGPos.y() - monoRhitGPos.y()) / sqrt(monoTsosGPEr.cyy() + monoRhitGPEr.cyy());
456  double pullGPZ_tr_mono = (monoTsosGPos.z() - monoRhitGPos.z()) / sqrt(monoTsosGPEr.czz() + monoRhitGPEr.czz());
457  //double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x());
458  //double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y());
459  //double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z());
460 
461  title.str("");
462  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
463  hPullGP_X_tr_mono[title.str()]->Fill(pullGPX_tr_mono);
464  title.str("");
465  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
466  hPullGP_Y_tr_mono[title.str()]->Fill(pullGPY_tr_mono);
467  title.str("");
468  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
469  hPullGP_Z_tr_mono[title.str()]->Fill(pullGPZ_tr_mono);
470 
471  double pullGPX_ts_mono = (monoTsosGPos.x() - monoShitGPos.x()) / sqrt(monoTsosGPEr.cxx());
472  double pullGPY_ts_mono = (monoTsosGPos.y() - monoShitGPos.y()) / sqrt(monoTsosGPEr.cyy());
473  double pullGPZ_ts_mono = (monoTsosGPos.z() - monoShitGPos.z()) / sqrt(monoTsosGPEr.czz());
474  //double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x());
475  //double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y());
476  //double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z());
477 
478  title.str("");
479  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
480  hPullGP_X_ts_mono[title.str()]->Fill(pullGPX_ts_mono);
481  title.str("");
482  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
483  hPullGP_Y_ts_mono[title.str()]->Fill(pullGPY_ts_mono);
484  title.str("");
485  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
486  hPullGP_Z_ts_mono[title.str()]->Fill(pullGPZ_ts_mono);
487 
488  double pullGMX_ts_mono = (monoTsosGMom.x() - monoShitGMom.x()) / sqrt(monoTsosGMEr.cxx());
489  double pullGMY_ts_mono = (monoTsosGMom.y() - monoShitGMom.y()) / sqrt(monoTsosGMEr.cyy());
490  double pullGMZ_ts_mono = (monoTsosGMom.z() - monoShitGMom.z()) / sqrt(monoTsosGMEr.czz());
491  //double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x());
492  //double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y());
493  //double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z());
494 
495  title.str("");
496  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
497  hPullGM_X_ts_mono[title.str()]->Fill(pullGMX_ts_mono);
498  title.str("");
499  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
500  hPullGM_Y_ts_mono[title.str()]->Fill(pullGMY_ts_mono);
501  title.str("");
502  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
503  hPullGM_Z_ts_mono[title.str()]->Fill(pullGMZ_ts_mono);
504 
505  //stereo
506  LogTrace("TestHits") << "STEREO HIT" << std::endl;
507  auto s = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->stereoHit();
508  CTTRHp tStereoHit = theBuilder->build(&s);
509  if (tStereoHit == nullptr)
510  continue;
511  vector<PSimHit> assStereoSimHits = hitAssociator.associateHit(*tStereoHit->hit());
512  if (assStereoSimHits.empty())
513  continue;
514  const PSimHit sStereoHit = *(assSimHits.begin());
515  const Surface* stereoSurf = &(tStereoHit->det()->surface());
516  if (stereoSurf == nullptr)
517  continue;
518  TSOS stereoState = thePropagator->propagate(lastState, *stereoSurf);
519  if (stereoState.isValid() == 0)
520  continue;
521 
522  LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
523  GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
524  LocalPoint stereoShitLPos = sStereoHit.localPosition();
525  GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
526 
527  GlobalVector stereoTsosGMom = stereoState.globalMomentum();
528  GlobalError stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
529  GlobalPoint stereoTsosGPos = stereoState.globalPosition();
530  GlobalError stereoTsosGPEr = stereoState.cartesianError().position();
531 
532  GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
533  GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
534 
535  double pullGPX_rs_stereo = (stereoRhitGPos.x() - stereoShitGPos.x()) / sqrt(stereoRhitGPEr.cxx());
536  double pullGPY_rs_stereo = (stereoRhitGPos.y() - stereoShitGPos.y()) / sqrt(stereoRhitGPEr.cyy());
537  double pullGPZ_rs_stereo = (stereoRhitGPos.z() - stereoShitGPos.z()) / sqrt(stereoRhitGPEr.czz());
538  //double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x());
539  //double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y());
540  //double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z());
541 
542  title.str("");
543  title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
544  hPullGP_X_rs_stereo[title.str()]->Fill(pullGPX_rs_stereo);
545  title.str("");
546  title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
547  hPullGP_Y_rs_stereo[title.str()]->Fill(pullGPY_rs_stereo);
548  title.str("");
549  title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
550  hPullGP_Z_rs_stereo[title.str()]->Fill(pullGPZ_rs_stereo);
551 
552  double pullGPX_tr_stereo =
553  (stereoTsosGPos.x() - stereoRhitGPos.x()) / sqrt(stereoTsosGPEr.cxx() + stereoRhitGPEr.cxx());
554  double pullGPY_tr_stereo =
555  (stereoTsosGPos.y() - stereoRhitGPos.y()) / sqrt(stereoTsosGPEr.cyy() + stereoRhitGPEr.cyy());
556  double pullGPZ_tr_stereo =
557  (stereoTsosGPos.z() - stereoRhitGPos.z()) / sqrt(stereoTsosGPEr.czz() + stereoRhitGPEr.czz());
558  //double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x());
559  //double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y());
560  //double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z());
561 
562  title.str("");
563  title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
564  hPullGP_X_tr_stereo[title.str()]->Fill(pullGPX_tr_stereo);
565  title.str("");
566  title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
567  hPullGP_Y_tr_stereo[title.str()]->Fill(pullGPY_tr_stereo);
568  title.str("");
569  title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
570  hPullGP_Z_tr_stereo[title.str()]->Fill(pullGPZ_tr_stereo);
571 
572  double pullGPX_ts_stereo = (stereoTsosGPos.x() - stereoShitGPos.x()) / sqrt(stereoTsosGPEr.cxx());
573  double pullGPY_ts_stereo = (stereoTsosGPos.y() - stereoShitGPos.y()) / sqrt(stereoTsosGPEr.cyy());
574  double pullGPZ_ts_stereo = (stereoTsosGPos.z() - stereoShitGPos.z()) / sqrt(stereoTsosGPEr.czz());
575  //double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x());
576  //double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y());
577  //double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z());
578 
579  title.str("");
580  title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
581  hPullGP_X_ts_stereo[title.str()]->Fill(pullGPX_ts_stereo);
582  title.str("");
583  title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
584  hPullGP_Y_ts_stereo[title.str()]->Fill(pullGPY_ts_stereo);
585  title.str("");
586  title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
587  hPullGP_Z_ts_stereo[title.str()]->Fill(pullGPZ_ts_stereo);
588 
589  double pullGMX_ts_stereo = (stereoTsosGMom.x() - stereoShitGMom.x()) / sqrt(stereoTsosGMEr.cxx());
590  double pullGMY_ts_stereo = (stereoTsosGMom.y() - stereoShitGMom.y()) / sqrt(stereoTsosGMEr.cyy());
591  double pullGMZ_ts_stereo = (stereoTsosGMom.z() - stereoShitGMom.z()) / sqrt(stereoTsosGMEr.czz());
592  //double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x());
593  //double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y());
594  //double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z());
595 
596  title.str("");
597  title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
598  hPullGM_X_ts_stereo[title.str()]->Fill(pullGMX_ts_stereo);
599  title.str("");
600  title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
601  hPullGM_Y_ts_stereo[title.str()]->Fill(pullGMY_ts_stereo);
602  title.str("");
603  title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
604  hPullGM_Z_ts_stereo[title.str()]->Fill(pullGMZ_ts_stereo);
605  }
606  lastState = updatedState;
607  }
608  LogTrace("TestHits") << "traj chi2=" << tchi2;
609  LogTrace("TestHits") << "track chi2=" << result[0].chiSquared();
610  }
611  LogTrace("TestHits") << "end of event" << std::endl;
612 }
613 // TSOS lastState = theTSOS;
614 // for (std::vector<TrajectoryMeasurement>::iterator tm=vtm.begin(); tm!=vtm.end();tm++){
615 
616 // TransientTrackingRecHit::ConstRecHitPointer rhit = tm->recHit();
617 // if ((rhit)->isValid()==0&&rhit->det()!=0) continue;
618 
619 // // //test hits
620 // // TSOS lastState, currentState, updatedState;
621 // // updatedState=theTSOS;
622 // // for (TransientTrackingRecHit::RecHitContainer::iterator rhit=hits.begin()+1;rhit!=hits.end();++rhit){
623 
624 // lastState=updatedState;
625 
626 // LogTrace("TestHits") << "new hit" << std::endl;
627 
628 // if ((*rhit)->isValid()==0) continue;
629 
630 // int subdetId = (*rhit)->det()->geographicalId().subdetId();
631 // int layerId = 0;
632 // DetId id = (*rhit)->det()->geographicalId();
633 // if (id.subdetId()==3) layerId = ((tTopo->tibLayer(id);
634 // if (id.subdetId()==5) layerId = ((tTopo->tobLayer(id);
635 // if (id.subdetId()==1) layerId = ((tTopo->pxbLayer(id);
636 // if (id.subdetId()==4) layerId = ((tTopo->tidWheel(id);
637 // if (id.subdetId()==6) layerId = ((tTopo->tecWheel(id);
638 // if (id.subdetId()==2) layerId = ((tTopo->pxfDisk(id);
639 // const Surface * surf = &( (*rhit)->det()->surface() );
640 // currentState=thePropagator->propagate(lastState,*surf);
641 // if (currentState.isValid()==0) continue;
642 // updatedState=theUpdator->update(currentState,**rhit);
643 
644 // double delta = 99999;
645 // LocalPoint rhitLP = rhit->localPosition();
646 
647 // std::vector<PSimHit> assSimHits = hitAssociator.associateHit(*(*rhit)->hit());
648 // if (assSimHits.size()==0) continue;
649 // PSimHit shit;
650 // for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
651 // if ((*m-rhitLP).mag()<delta) {
652 // shit=*m;
653 // delta = (*m-rhitLP).mag();
654 // }
655 // }
656 
657 // LocalVector shitLMom = shit.momentumAtEntry();
658 // GlobalVector shitGMom = surf->toGlobal(shitLMom);
659 // LocalPoint shitLPos = shit.localPosition();
660 // GlobalPoint shitGPos = surf->toGlobal(shitLPos);
661 
662 // GlobalVector tsosGMom = currentState.globalMomentum();
663 // GlobalError tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
664 // GlobalPoint tsosGPos = currentState.globalPosition();
665 // GlobalError tsosGPEr = currentState.cartesianError().position();
666 
667 // GlobalPoint rhitGPos = (*rhit)->globalPosition();
668 // GlobalError rhitGPEr = (*rhit)->globalPositionError();
669 
670 // double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/sqrt(rhitGPEr.cxx());
671 // double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/sqrt(rhitGPEr.cyy());
672 // double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/sqrt(rhitGPEr.czz());
673 // // double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/*/sqrt(rhitGPEr.cxx())*/;
674 // // double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/*/sqrt(rhitGPEr.cyy())*/;
675 // // double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/*/sqrt(rhitGPEr.czz())*/;
676 
677 // //plot chi2 increment
678 // MeasurementExtractor me(currentState);
679 // double chi2increment = computeChi2Increment(me,*rhit);
680 // LogTrace("TestHits") << "chi2increment=" << chi2increment << std::endl;
681 // title.str("");
682 // title << "Chi2Increment_" << subdetId << "-" << layerId;
683 // hChi2Increment[title.str()]->Fill( chi2increment );
684 // hTotChi2Increment->Fill( chi2increment );
685 // hChi2_vs_Process->Fill( chi2increment, shit.processType() );
686 // if (dynamic_cast<const SiPixelRecHit*>((*rhit)->hit()))
687 // hChi2_vs_clsize->Fill( chi2increment, ((const SiPixelRecHit*)(*rhit)->hit())->cluster()->size() );
688 // if (dynamic_cast<const SiStripRecHit2D*>((*rhit)->hit()))
689 // hChi2_vs_clsize->Fill( chi2increment, ((const SiStripRecHit2D*)(*rhit)->hit())->cluster()->amplitudes().size() );
690 
691 // LogTrace("TestHits") << "rs" << std::endl;
692 
693 // title.str("");
694 // title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
695 // hPullGP_X_rs[title.str()]->Fill( pullGPX_rs );
696 // title.str("");
697 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
698 // hPullGP_Y_rs[title.str()]->Fill( pullGPY_rs );
699 // title.str("");
700 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
701 // hPullGP_Z_rs[title.str()]->Fill( pullGPZ_rs );
702 
703 // double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx());
704 // double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy());
705 // double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/sqrt(tsosGPEr.czz()+rhitGPEr.czz());
706 // // double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/*/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx())*/;
707 // // double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/*/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy())*/;
708 // // double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/*/sqrt(tsosGPEr.czz()+rhitGPEr.czz())*/;
709 
710 // LogTrace("TestHits") << "tr" << std::endl;
711 
712 // title.str("");
713 // title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
714 // hPullGP_X_tr[title.str()]->Fill( pullGPX_tr );
715 // title.str("");
716 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
717 // hPullGP_Y_tr[title.str()]->Fill( pullGPY_tr );
718 // title.str("");
719 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
720 // hPullGP_Z_tr[title.str()]->Fill( pullGPZ_tr );
721 
722 // double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/sqrt(tsosGPEr.cxx());
723 // double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/sqrt(tsosGPEr.cyy());
724 // double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/sqrt(tsosGPEr.czz());
725 // // double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/*/sqrt(tsosGPEr.cxx())*/;
726 // // double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/*/sqrt(tsosGPEr.cyy())*/;
727 // // double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/*/sqrt(tsosGPEr.czz())*/;
728 
729 // LogTrace("TestHits") << "ts1" << std::endl;
730 
731 // title.str("");
732 // title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
733 // hPullGP_X_ts[title.str()]->Fill( pullGPX_ts );
734 // title.str("");
735 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
736 // hPullGP_Y_ts[title.str()]->Fill( pullGPY_ts );
737 // title.str("");
738 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
739 // hPullGP_Z_ts[title.str()]->Fill( pullGPZ_ts );
740 
741 // double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/sqrt(tsosGMEr.cxx());
742 // double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/sqrt(tsosGMEr.cyy());
743 // double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/sqrt(tsosGMEr.czz());
744 // // double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/*/sqrt(tsosGMEr.cxx())*/;
745 // // double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/*/sqrt(tsosGMEr.cyy())*/;
746 // // double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/*/sqrt(tsosGMEr.czz())*/;
747 
748 // LogTrace("TestHits") << "ts2" << std::endl;
749 
750 // title.str("");
751 // title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
752 // hPullGM_X_ts[title.str()]->Fill( pullGMX_ts );
753 // title.str("");
754 // title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
755 // hPullGM_Y_ts[title.str()]->Fill( pullGMY_ts );
756 // title.str("");
757 // title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
758 // hPullGM_Z_ts[title.str()]->Fill( pullGMZ_ts );
759 
760 // if (dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())) {
761 // //mono
762 // LogTrace("TestHits") << "MONO HIT" << std::endl;
763 // CTTRHp tMonoHit =
764 // theBuilder->build(dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())->monoHit());
765 // if (tMonoHit==0) continue;
766 // vector<PSimHit> assMonoSimHits = hitAssociator.associateHit(*tMonoHit->hit());
767 // if (assMonoSimHits.size()==0) continue;
768 // const PSimHit sMonoHit = *(assSimHits.begin());
769 // const Surface * monoSurf = &( tMonoHit->det()->surface() );
770 // if (monoSurf==0) continue;
771 // TSOS monoState = thePropagator->propagate(lastState,*monoSurf);
772 // if (monoState.isValid()==0) continue;
773 
774 // LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
775 // GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
776 // LocalPoint monoShitLPos = sMonoHit.localPosition();
777 // GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
778 
779 // GlobalVector monoTsosGMom = monoState.globalMomentum();
780 // GlobalError monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
781 // GlobalPoint monoTsosGPos = monoState.globalPosition();
782 // GlobalError monoTsosGPEr = monoState.cartesianError().position();
783 
784 // GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
785 // GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
786 
787 // double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/sqrt(monoRhitGPEr.cxx());
788 // double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/sqrt(monoRhitGPEr.cyy());
789 // double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/sqrt(monoRhitGPEr.czz());
790 // // double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/*/sqrt(monoRhitGPEr.cxx())*/;
791 // // double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/*/sqrt(monoRhitGPEr.cyy())*/;
792 // // double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/*/sqrt(monoRhitGPEr.czz())*/;
793 
794 // title.str("");
795 // title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
796 // hPullGP_X_rs_mono[title.str()]->Fill( pullGPX_rs_mono );
797 // title.str("");
798 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
799 // hPullGP_Y_rs_mono[title.str()]->Fill( pullGPY_rs_mono );
800 // title.str("");
801 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
802 // hPullGP_Z_rs_mono[title.str()]->Fill( pullGPZ_rs_mono );
803 
804 // double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx());
805 // double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy());
806 // double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz());
807 // // double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/*/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx())*/;
808 // // double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/*/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy())*/;
809 // // double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/*/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz())*/;
810 
811 // title.str("");
812 // title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
813 // hPullGP_X_tr_mono[title.str()]->Fill( pullGPX_tr_mono );
814 // title.str("");
815 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
816 // hPullGP_Y_tr_mono[title.str()]->Fill( pullGPY_tr_mono );
817 // title.str("");
818 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
819 // hPullGP_Z_tr_mono[title.str()]->Fill( pullGPZ_tr_mono );
820 
821 // double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/sqrt(monoTsosGPEr.cxx());
822 // double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/sqrt(monoTsosGPEr.cyy());
823 // double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/sqrt(monoTsosGPEr.czz());
824 // // double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/*/sqrt(monoTsosGPEr.cxx())*/;
825 // // double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/*/sqrt(monoTsosGPEr.cyy())*/;
826 // // double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/*/sqrt(monoTsosGPEr.czz())*/;
827 
828 // title.str("");
829 // title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
830 // hPullGP_X_ts_mono[title.str()]->Fill( pullGPX_ts_mono );
831 // title.str("");
832 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
833 // hPullGP_Y_ts_mono[title.str()]->Fill( pullGPY_ts_mono );
834 // title.str("");
835 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
836 // hPullGP_Z_ts_mono[title.str()]->Fill( pullGPZ_ts_mono );
837 
838 // double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/sqrt(monoTsosGMEr.cxx());
839 // double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/sqrt(monoTsosGMEr.cyy());
840 // double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/sqrt(monoTsosGMEr.czz());
841 // // double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/*/sqrt(monoTsosGMEr.cxx())*/;
842 // // double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/*/sqrt(monoTsosGMEr.cyy())*/;
843 // // double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/*/sqrt(monoTsosGMEr.czz())*/;
844 
845 // title.str("");
846 // title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
847 // hPullGM_X_ts_mono[title.str()]->Fill( pullGMX_ts_mono );
848 // title.str("");
849 // title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
850 // hPullGM_Y_ts_mono[title.str()]->Fill( pullGMY_ts_mono );
851 // title.str("");
852 // title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
853 // hPullGM_Z_ts_mono[title.str()]->Fill( pullGMZ_ts_mono );
854 
855 // //stereo
856 // LogTrace("TestHits") << "STEREO HIT" << std::endl;
857 // CTTRHp tStereoHit =
858 // theBuilder->build(dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())->stereoHit());
859 // if (tStereoHit==0) continue;
860 // vector<PSimHit> assStereoSimHits = hitAssociator.associateHit(*tStereoHit->hit());
861 // if (assStereoSimHits.size()==0) continue;
862 // const PSimHit sStereoHit = *(assSimHits.begin());
863 // const Surface * stereoSurf = &( tStereoHit->det()->surface() );
864 // if (stereoSurf==0) continue;
865 // TSOS stereoState = thePropagator->propagate(lastState,*stereoSurf);
866 // if (stereoState.isValid()==0) continue;
867 
868 // LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
869 // GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
870 // LocalPoint stereoShitLPos = sStereoHit.localPosition();
871 // GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
872 
873 // GlobalVector stereoTsosGMom = stereoState.globalMomentum();
874 // GlobalError stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
875 // GlobalPoint stereoTsosGPos = stereoState.globalPosition();
876 // GlobalError stereoTsosGPEr = stereoState.cartesianError().position();
877 
878 // GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
879 // GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
880 
881 // double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/sqrt(stereoRhitGPEr.cxx());
882 // double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/sqrt(stereoRhitGPEr.cyy());
883 // double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/sqrt(stereoRhitGPEr.czz());
884 // // double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/*/sqrt(stereoRhitGPEr.cxx())*/;
885 // // double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/*/sqrt(stereoRhitGPEr.cyy())*/;
886 // // double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/*/sqrt(stereoRhitGPEr.czz())*/;
887 
888 // title.str("");
889 // title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
890 // hPullGP_X_rs_stereo[title.str()]->Fill( pullGPX_rs_stereo );
891 // title.str("");
892 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
893 // hPullGP_Y_rs_stereo[title.str()]->Fill( pullGPY_rs_stereo );
894 // title.str("");
895 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
896 // hPullGP_Z_rs_stereo[title.str()]->Fill( pullGPZ_rs_stereo );
897 
898 // double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx());
899 // double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy());
900 // double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz());
901 // // double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/*/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx())*/;
902 // // double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/*/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy())*/;
903 // // double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/*/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz())*/;
904 
905 // title.str("");
906 // title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
907 // hPullGP_X_tr_stereo[title.str()]->Fill( pullGPX_tr_stereo );
908 // title.str("");
909 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
910 // hPullGP_Y_tr_stereo[title.str()]->Fill( pullGPY_tr_stereo );
911 // title.str("");
912 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
913 // hPullGP_Z_tr_stereo[title.str()]->Fill( pullGPZ_tr_stereo );
914 
915 // double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/sqrt(stereoTsosGPEr.cxx());
916 // double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/sqrt(stereoTsosGPEr.cyy());
917 // double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/sqrt(stereoTsosGPEr.czz());
918 // // double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/*/sqrt(stereoTsosGPEr.cxx())*/;
919 // // double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/*/sqrt(stereoTsosGPEr.cyy())*/;
920 // // double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/*/sqrt(stereoTsosGPEr.czz())*/;
921 
922 // title.str("");
923 // title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
924 // hPullGP_X_ts_stereo[title.str()]->Fill( pullGPX_ts_stereo );
925 // title.str("");
926 // title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
927 // hPullGP_Y_ts_stereo[title.str()]->Fill( pullGPY_ts_stereo );
928 // title.str("");
929 // title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
930 // hPullGP_Z_ts_stereo[title.str()]->Fill( pullGPZ_ts_stereo );
931 
932 // double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/sqrt(stereoTsosGMEr.cxx());
933 // double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/sqrt(stereoTsosGMEr.cyy());
934 // double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/sqrt(stereoTsosGMEr.czz());
935 // // double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/*/sqrt(stereoTsosGMEr.cxx())*/;
936 // // double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/*/sqrt(stereoTsosGMEr.cyy())*/;
937 // // double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/*/sqrt(stereoTsosGMEr.czz())*/;
938 
939 // title.str("");
940 // title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
941 // hPullGM_X_ts_stereo[title.str()]->Fill( pullGMX_ts_stereo );
942 // title.str("");
943 // title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
944 // hPullGM_Y_ts_stereo[title.str()]->Fill( pullGMY_ts_stereo );
945 // title.str("");
946 // title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
947 // hPullGM_Z_ts_stereo[title.str()]->Fill( pullGMZ_ts_stereo );
948 // }
949 // }
950 // }
951 // LogTrace("TestHits") << "end of event" << std::endl;
952 // }
953 
955  //file->Write();
956  TDirectory* chi2i = file->mkdir("Chi2_Increment");
957 
958  TDirectory* gp_ts = file->mkdir("GP_TSOS-SimHit");
959  TDirectory* gm_ts = file->mkdir("GM_TSOS-SimHit");
960  TDirectory* gp_tr = file->mkdir("GP_TSOS-RecHit");
961  TDirectory* gp_rs = file->mkdir("GP_RecHit-SimHit");
962 
963  TDirectory* gp_tsx = gp_ts->mkdir("X");
964  TDirectory* gp_tsy = gp_ts->mkdir("Y");
965  TDirectory* gp_tsz = gp_ts->mkdir("Z");
966  TDirectory* gm_tsx = gm_ts->mkdir("X");
967  TDirectory* gm_tsy = gm_ts->mkdir("Y");
968  TDirectory* gm_tsz = gm_ts->mkdir("Z");
969  TDirectory* gp_trx = gp_tr->mkdir("X");
970  TDirectory* gp_try = gp_tr->mkdir("Y");
971  TDirectory* gp_trz = gp_tr->mkdir("Z");
972  TDirectory* gp_rsx = gp_rs->mkdir("X");
973  TDirectory* gp_rsy = gp_rs->mkdir("Y");
974  TDirectory* gp_rsz = gp_rs->mkdir("Z");
975 
976  TDirectory* gp_tsx_mono = gp_ts->mkdir("MONOX");
977  TDirectory* gp_tsy_mono = gp_ts->mkdir("MONOY");
978  TDirectory* gp_tsz_mono = gp_ts->mkdir("MONOZ");
979  TDirectory* gm_tsx_mono = gm_ts->mkdir("MONOX");
980  TDirectory* gm_tsy_mono = gm_ts->mkdir("MONOY");
981  TDirectory* gm_tsz_mono = gm_ts->mkdir("MONOZ");
982  TDirectory* gp_trx_mono = gp_tr->mkdir("MONOX");
983  TDirectory* gp_try_mono = gp_tr->mkdir("MONOY");
984  TDirectory* gp_trz_mono = gp_tr->mkdir("MONOZ");
985  TDirectory* gp_rsx_mono = gp_rs->mkdir("MONOX");
986  TDirectory* gp_rsy_mono = gp_rs->mkdir("MONOY");
987  TDirectory* gp_rsz_mono = gp_rs->mkdir("MONOZ");
988 
989  TDirectory* gp_tsx_stereo = gp_ts->mkdir("STEREOX");
990  TDirectory* gp_tsy_stereo = gp_ts->mkdir("STEREOY");
991  TDirectory* gp_tsz_stereo = gp_ts->mkdir("STEREOZ");
992  TDirectory* gm_tsx_stereo = gm_ts->mkdir("STEREOX");
993  TDirectory* gm_tsy_stereo = gm_ts->mkdir("STEREOY");
994  TDirectory* gm_tsz_stereo = gm_ts->mkdir("STEREOZ");
995  TDirectory* gp_trx_stereo = gp_tr->mkdir("STEREOX");
996  TDirectory* gp_try_stereo = gp_tr->mkdir("STEREOY");
997  TDirectory* gp_trz_stereo = gp_tr->mkdir("STEREOZ");
998  TDirectory* gp_rsx_stereo = gp_rs->mkdir("STEREOX");
999  TDirectory* gp_rsy_stereo = gp_rs->mkdir("STEREOY");
1000  TDirectory* gp_rsz_stereo = gp_rs->mkdir("STEREOZ");
1001 
1002  chi2i->cd();
1003  hTotChi2Increment->Write();
1004  hProcess_vs_Chi2->Write();
1005  hClsize_vs_Chi2->Write();
1006  for (int i = 0; i != 6; i++)
1007  for (int j = 0; j != 9; j++) {
1008  if (i == 0 && j > 2)
1009  break;
1010  if (i == 1 && j > 1)
1011  break;
1012  if (i == 2 && j > 3)
1013  break;
1014  if (i == 3 && j > 2)
1015  break;
1016  if (i == 4 && j > 5)
1017  break;
1018  if (i == 5 && j > 8)
1019  break;
1020  chi2i->cd();
1021  title.str("");
1022  title << "Chi2Increment_" << i + 1 << "-" << j + 1;
1023  hChi2Increment[title.str()]->Write();
1024 
1025  gp_ts->cd();
1026  gp_tsx->cd();
1027  title.str("");
1028  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts";
1029  hPullGP_X_ts[title.str()]->Write();
1030  gp_tsy->cd();
1031  title.str("");
1032  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts";
1033  hPullGP_Y_ts[title.str()]->Write();
1034  gp_tsz->cd();
1035  title.str("");
1036  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts";
1037  hPullGP_Z_ts[title.str()]->Write();
1038 
1039  gm_ts->cd();
1040  gm_tsx->cd();
1041  title.str("");
1042  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts";
1043  hPullGM_X_ts[title.str()]->Write();
1044  gm_tsy->cd();
1045  title.str("");
1046  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts";
1047  hPullGM_Y_ts[title.str()]->Write();
1048  gm_tsz->cd();
1049  title.str("");
1050  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts";
1051  hPullGM_Z_ts[title.str()]->Write();
1052 
1053  gp_tr->cd();
1054  gp_trx->cd();
1055  title.str("");
1056  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr";
1057  hPullGP_X_tr[title.str()]->Write();
1058  gp_try->cd();
1059  title.str("");
1060  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr";
1061  hPullGP_Y_tr[title.str()]->Write();
1062  gp_trz->cd();
1063  title.str("");
1064  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr";
1065  hPullGP_Z_tr[title.str()]->Write();
1066 
1067  gp_rs->cd();
1068  gp_rsx->cd();
1069  title.str("");
1070  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs";
1071  hPullGP_X_rs[title.str()]->Write();
1072  gp_rsy->cd();
1073  title.str("");
1074  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs";
1075  hPullGP_Y_rs[title.str()]->Write();
1076  gp_rsz->cd();
1077  title.str("");
1078  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs";
1079  hPullGP_Z_rs[title.str()]->Write();
1080 
1081  if (((i == 2 || i == 4) && (j == 0 || j == 1)) || (i == 3 || i == 5)) {
1082  //mono
1083  gp_ts->cd();
1084  gp_tsx_mono->cd();
1085  title.str("");
1086  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
1087  hPullGP_X_ts_mono[title.str()]->Write();
1088  gp_tsy_mono->cd();
1089  title.str("");
1090  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
1091  hPullGP_Y_ts_mono[title.str()]->Write();
1092  gp_tsz_mono->cd();
1093  title.str("");
1094  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
1095  hPullGP_Z_ts_mono[title.str()]->Write();
1096 
1097  gm_ts->cd();
1098  gm_tsx_mono->cd();
1099  title.str("");
1100  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
1101  hPullGM_X_ts_mono[title.str()]->Write();
1102  gm_tsy_mono->cd();
1103  title.str("");
1104  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
1105  hPullGM_Y_ts_mono[title.str()]->Write();
1106  gm_tsz_mono->cd();
1107  title.str("");
1108  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
1109  hPullGM_Z_ts_mono[title.str()]->Write();
1110 
1111  gp_tr->cd();
1112  gp_trx_mono->cd();
1113  title.str("");
1114  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_mono";
1115  hPullGP_X_tr_mono[title.str()]->Write();
1116  gp_try_mono->cd();
1117  title.str("");
1118  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_mono";
1119  hPullGP_Y_tr_mono[title.str()]->Write();
1120  gp_trz_mono->cd();
1121  title.str("");
1122  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_mono";
1123  hPullGP_Z_tr_mono[title.str()]->Write();
1124 
1125  gp_rs->cd();
1126  gp_rsx_mono->cd();
1127  title.str("");
1128  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_mono";
1129  hPullGP_X_rs_mono[title.str()]->Write();
1130  gp_rsy_mono->cd();
1131  title.str("");
1132  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_mono";
1133  hPullGP_Y_rs_mono[title.str()]->Write();
1134  gp_rsz_mono->cd();
1135  title.str("");
1136  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_mono";
1137  hPullGP_Z_rs_mono[title.str()]->Write();
1138 
1139  //stereo
1140  gp_ts->cd();
1141  gp_tsx_stereo->cd();
1142  title.str("");
1143  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1144  hPullGP_X_ts_stereo[title.str()]->Write();
1145  gp_tsy_stereo->cd();
1146  title.str("");
1147  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1148  hPullGP_Y_ts_stereo[title.str()]->Write();
1149  gp_tsz_stereo->cd();
1150  title.str("");
1151  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1152  hPullGP_Z_ts_stereo[title.str()]->Write();
1153 
1154  gm_ts->cd();
1155  gm_tsx_stereo->cd();
1156  title.str("");
1157  title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1158  hPullGM_X_ts_stereo[title.str()]->Write();
1159  gm_tsy_stereo->cd();
1160  title.str("");
1161  title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1162  hPullGM_Y_ts_stereo[title.str()]->Write();
1163  gm_tsz_stereo->cd();
1164  title.str("");
1165  title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
1166  hPullGM_Z_ts_stereo[title.str()]->Write();
1167 
1168  gp_tr->cd();
1169  gp_trx_stereo->cd();
1170  title.str("");
1171  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_stereo";
1172  hPullGP_X_tr_stereo[title.str()]->Write();
1173  gp_try_stereo->cd();
1174  title.str("");
1175  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_stereo";
1176  hPullGP_Y_tr_stereo[title.str()]->Write();
1177  gp_trz_stereo->cd();
1178  title.str("");
1179  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_stereo";
1180  hPullGP_Z_tr_stereo[title.str()]->Write();
1181 
1182  gp_rs->cd();
1183  gp_rsx_stereo->cd();
1184  title.str("");
1185  title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_stereo";
1186  hPullGP_X_rs_stereo[title.str()]->Write();
1187  gp_rsy_stereo->cd();
1188  title.str("");
1189  title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_stereo";
1190  hPullGP_Y_rs_stereo[title.str()]->Write();
1191  gp_rsz_stereo->cd();
1192  title.str("");
1193  title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_stereo";
1194  hPullGP_Z_rs_stereo[title.str()]->Write();
1195  }
1196  }
1197 
1198  file->Close();
1199 }
1200 
1201 //needed by to do the residual for matched hits
1202 //taken from SiStripTrackingRecHitsValid.cc
1203 std::pair<LocalPoint, LocalVector> TestHits::projectHit(const PSimHit& hit,
1204  const StripGeomDetUnit* stripDet,
1205  const BoundPlane& plane) {
1206  const StripTopology& topol = stripDet->specificTopology();
1207  GlobalPoint globalpos = stripDet->surface().toGlobal(hit.localPosition());
1208  LocalPoint localHit = plane.toLocal(globalpos);
1209  //track direction
1210  LocalVector locdir = hit.localDirection();
1211  //rotate track in new frame
1212 
1213  GlobalVector globaldir = stripDet->surface().toGlobal(locdir);
1214  LocalVector dir = plane.toLocal(globaldir);
1215  float scale = -localHit.z() / dir.z();
1216 
1217  LocalPoint projectedPos = localHit + scale * dir;
1218 
1219  float selfAngle = topol.stripAngle(topol.strip(hit.localPosition()));
1220 
1221  LocalVector stripDir(sin(selfAngle), cos(selfAngle), 0); // vector along strip in hit frame
1222 
1223  LocalVector localStripDir(plane.toLocal(stripDet->surface().toGlobal(stripDir)));
1224 
1225  return std::pair<LocalPoint, LocalVector>(projectedPos, localStripDir);
1226 }
1227 
1230 
edm::Handle< TrackCandidateCollection > theTCCollection
Definition: TestHits.h:76
Log< level::Info, true > LogVerbatim
std::map< std::string, TH1F * > hPullGP_Y_tr_mono
Definition: TestHits.h:114
std::map< std::string, TH1F * > hPullGM_Z_ts_stereo
Definition: TestHits.h:122
std::map< std::string, TH1F * > hPullGM_X_ts_mono
Definition: TestHits.h:107
std::map< std::string, TH1F * > hPullGP_X_tr
Definition: TestHits.h:97
std::pair< LocalPoint, LocalVector > projectHit(const PSimHit &, const StripGeomDetUnit *, const BoundPlane &)
Definition: TestHits.cc:1203
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
std::map< std::string, TH1F * > hPullGP_X_ts
Definition: TestHits.h:88
std::map< std::string, TH1F * > hPullGM_Z_ts
Definition: TestHits.h:93
std::map< std::string, TH1F * > hPullGP_X_rs_mono
Definition: TestHits.h:110
TH1F * hTotChi2Increment
Definition: TestHits.h:101
std::map< std::string, TH1F * > hPullGP_Z_tr
Definition: TestHits.h:99
std::map< std::string, TH1F * > hPullGM_X_ts_stereo
Definition: TestHits.h:120
const bool isValid(const Frame &aFrame, const FrameQuality &aQuality, const uint16_t aExpectedPos)
virtual float stripAngle(float strip) const =0
std::map< std::string, TH1F * > hPullGP_Z_tr_mono
Definition: TestHits.h:115
edm::ESHandle< TransientTrackingRecHitBuilder > theBuilder
Definition: TestHits.h:74
std::map< std::string, TH1F * > hPullGM_Z_ts_mono
Definition: TestHits.h:109
T z() const
Definition: PV3DBase.h:61
std::map< std::string, TH1F * > hPullGP_Y_rs_stereo
Definition: TestHits.h:124
std::string builderName
Definition: TestHits.h:67
double mineta
Definition: TestHits.h:64
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
edm::ESGetToken< TrajectoryFitter, TrajectoryFitter::Record > fitToken
Definition: TestHits.h:82
std::map< std::string, TH1F * > hPullGP_Z_rs
Definition: TestHits.h:96
std::map< std::string, TH1F * > hPullGP_Z_rs_stereo
Definition: TestHits.h:125
std::vector< ConstRecHitPointer > RecHitContainer
void beginRun(edm::Run const &run, const edm::EventSetup &) override
Definition: TestHits.cc:45
edm::ESGetToken< MagneticField, IdealMagneticFieldRecord > theMFToken
Definition: TestHits.h:79
std::map< std::string, TH1F * > hPullGP_X_tr_mono
Definition: TestHits.h:113
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken
Definition: TestHits.h:84
TrackerHitAssociator::Config trackerHitAssociatorConfig_
Definition: TestHits.h:62
std::map< std::string, TH1F * > hPullGP_Y_rs_mono
Definition: TestHits.h:111
TrajectoryStateOnSurface propagate(STA const &state, SUR const &surface) const
Definition: Propagator.h:50
constexpr Detector det() const
get the detector field from this detid
Definition: DetId.h:46
edm::ESHandle< Propagator > thePropagator
Definition: TestHits.h:73
edm::ESGetToken< Propagator, TrackingComponentsRecord > thePropagatorToken
Definition: TestHits.h:80
const CartesianTrajectoryError cartesianError() const
unsigned int layer(const DetId &id) const
#define LogTrace(id)
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
virtual float strip(const LocalPoint &) const =0
std::map< std::string, TH1F * > hPullGP_Y_tr
Definition: TestHits.h:98
edm::ESGetToken< TransientTrackingRecHitBuilder, TransientRecHitRecord > theBuilderToken
Definition: TestHits.h:81
std::map< std::string, TH1F * > hPullGP_X_rs
Definition: TestHits.h:94
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
std::map< std::string, TH1F * > hPullGP_X_rs_stereo
Definition: TestHits.h:123
int iEvent
Definition: GenABIO.cc:224
T const * product() const
Definition: ESHandle.h:86
std::map< std::string, TH1F * > hPullGP_X_tr_stereo
Definition: TestHits.h:126
GlobalPoint globalPosition() const
std::map< std::string, TH1F * > hPullGP_X_ts_stereo
Definition: TestHits.h:117
std::string propagatorName
Definition: TestHits.h:66
std::shared_ptr< TrackingRecHit const > ConstRecHitPointer
const GlobalError position() const
Position error submatrix.
std::map< std::string, TH1F * > hPullGP_Y_ts
Definition: TestHits.h:89
T sqrt(T t)
Definition: SSEVec.h:19
virtual RecHitPointer build(const TrackingRecHit *p) const =0
build a tracking rechit from an existing rechit
unsigned short processType() const
Definition: PSimHit.h:131
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
PTrajectoryStateOnDet const & trajectoryStateOnDet() const
std::map< std::string, TH1F * > hPullGP_Y_ts_mono
Definition: TestHits.h:105
virtual const StripTopology & specificTopology() const
Returns a reference to the strip proxy topology.
TransientTrackingRecHit::ConstRecHitPointer CTTRHp
Definition: CkfDebugger.h:44
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
std::map< std::string, TH1F * > hPullGP_Z_ts_mono
Definition: TestHits.h:106
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:130
const TrackerGeomDet * idToDet(DetId) const override
edm::ESHandle< TrackerGeometry > theG
Definition: TestHits.h:71
std::map< std::string, TH1F * > hPullGP_Z_ts
Definition: TestHits.h:90
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:79
Definition: DetId.h:17
TrajectoryStateOnSurface TSOS
Definition: TestHits.cc:19
std::map< std::string, TH1F * > hPullGM_Y_ts_stereo
Definition: TestHits.h:121
std::map< std::string, TH1F * > hPullGP_Y_tr_stereo
Definition: TestHits.h:127
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
const PositionType & position() const
TrajectoryStateOnSurface transientState(const PTrajectoryStateOnDet &ts, const Surface *surface, const MagneticField *field)
edm::EDGetTokenT< TrackCandidateCollection > theTCCollectionToken
Definition: TestHits.h:83
~TestHits() override
Definition: TestHits.cc:43
TH2F * hProcess_vs_Chi2
Definition: TestHits.h:102
TransientTrackingRecHit::ConstRecHitPointer CTTRHp
Definition: TestHits.cc:20
Local3DPoint localPosition() const
Definition: PSimHit.h:54
std::map< std::string, TH1F * > hPullGP_Z_rs_mono
Definition: TestHits.h:112
std::map< std::string, TH1F * > hPullGM_Y_ts
Definition: TestHits.h:92
GlobalVector globalMomentum() const
edm::ESHandle< MagneticField > theMF
Definition: TestHits.h:72
std::map< std::string, TH1F * > hChi2Increment
Definition: TestHits.h:100
TH2F * hClsize_vs_Chi2
Definition: TestHits.h:102
std::map< std::string, TH1F * > hPullGM_X_ts
Definition: TestHits.h:91
std::map< std::string, TH1F * > hPullGP_Y_rs
Definition: TestHits.h:95
HLT enums.
TFile * file
Definition: TestHits.h:86
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > theGToken
Definition: TestHits.h:78
std::map< std::string, TH1F * > hPullGM_Y_ts_mono
Definition: TestHits.h:108
void analyze(const edm::Event &, const edm::EventSetup &) override
Definition: TestHits.cc:199
const AlgebraicSymMatrix66 & matrix() const
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
TrajectorySeed const & seed() const
std::string fname
Definition: TestHits.h:69
std::map< std::string, TH1F * > hPullGP_Y_ts_stereo
Definition: TestHits.h:118
std::string srcName
Definition: TestHits.h:68
std::map< std::string, TH1F * > hPullGP_Z_ts_stereo
Definition: TestHits.h:119
TestHits(const edm::ParameterSet &)
Definition: TestHits.cc:24
edm::Range< RecHitContainer::const_iterator > recHits() const
LocalVector momentumAtEntry() const
The momentum of the track that produced the hit, at entry point.
Definition: PSimHit.h:57
std::map< std::string, TH1F * > hPullGP_X_ts_mono
Definition: TestHits.h:104
Definition: Run.h:45
std::stringstream title
Definition: TestHits.h:87
std::map< std::string, TH1F * > hPullGP_Z_tr_stereo
Definition: TestHits.h:128
Our base class.
Definition: SiPixelRecHit.h:23
double maxeta
Definition: TestHits.h:64
void endJob() override
Definition: TestHits.cc:954