File indexing completed on 2024-04-06 12:28:04
0001 #include "TestSmoothHits.h"
0002
0003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0004 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0005 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
0006 #include "TrackingTools/PatternTools/interface/TransverseImpactPointExtrapolator.h"
0007 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0008 #include "DataFormats/TrackerRecHit2D/interface/TrackingRecHitLessFromGlobalPosition.h"
0009 #include "TrackingTools/PatternTools/interface/TSCPBuilderNoMaterial.h"
0010 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
0011 #include <TDirectory.h>
0012 #include "Geometry/CommonDetUnit/interface/GluedGeomDet.h"
0013
0014 #include "TrackingTools/TrackFitters/interface/TrajectoryFitter.h"
0015 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0016 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0017
0018 typedef TrajectoryStateOnSurface TSOS;
0019 typedef TransientTrackingRecHit::ConstRecHitPointer CTTRHp;
0020 using namespace std;
0021 using namespace edm;
0022
0023 TestSmoothHits::TestSmoothHits(const edm::ParameterSet& iConfig) : trackerHitAssociatorConfig_(consumesCollector()) {
0024 LogTrace("TestSmoothHits") << iConfig << std::endl;
0025 propagatorName = iConfig.getParameter<std::string>("Propagator");
0026 builderName = iConfig.getParameter<std::string>("TTRHBuilder");
0027 srcName = iConfig.getParameter<std::string>("src");
0028 fname = iConfig.getParameter<std::string>("Fitter");
0029 sname = iConfig.getParameter<std::string>("Smoother");
0030 mineta = iConfig.getParameter<double>("mineta");
0031 maxeta = iConfig.getParameter<double>("maxeta");
0032
0033 theGToken = esConsumes<edm::Transition::BeginRun>();
0034 theMFToken = esConsumes<edm::Transition::BeginRun>();
0035 thePropagatorToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", propagatorName));
0036 theBuilderToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", builderName));
0037 fitToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", fname));
0038 smoothToken = esConsumes<edm::Transition::BeginRun>(edm::ESInputTag("", sname));
0039 tTopoToken = esConsumes();
0040
0041 theTCCollectionToken = consumes(edm::InputTag(srcName));
0042 }
0043
0044 TestSmoothHits::~TestSmoothHits() {}
0045
0046 void TestSmoothHits::beginRun(edm::Run const& run, const edm::EventSetup& iSetup) {
0047 theG = iSetup.getHandle(theGToken);
0048 theMF = iSetup.getHandle(theMFToken);
0049 thePropagator = iSetup.getHandle(thePropagatorToken);
0050 theBuilder = iSetup.getHandle(theBuilderToken);
0051 fit = iSetup.getHandle(fitToken);
0052 smooth = iSetup.getHandle(smoothToken);
0053
0054 file = new TFile("testSmoothHits.root", "recreate");
0055 for (int i = 0; i != 6; i++)
0056 for (int j = 0; j != 9; j++) {
0057 if (i == 0 && j > 2)
0058 break;
0059 if (i == 1 && j > 1)
0060 break;
0061 if (i == 2 && j > 3)
0062 break;
0063 if (i == 3 && j > 2)
0064 break;
0065 if (i == 4 && j > 5)
0066 break;
0067 if (i == 5 && j > 8)
0068 break;
0069 title.str("");
0070 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts";
0071 hPullGP_X_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0072 title.str("");
0073 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts";
0074 hPullGP_Y_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0075 title.str("");
0076 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts";
0077 hPullGP_Z_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0078 title.str("");
0079 title << "Chi2Increment_" << i + 1 << "-" << j + 1;
0080 hChi2Increment[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, 0, 100);
0081
0082 title.str("");
0083 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts";
0084 hPullGM_X_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0085 title.str("");
0086 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts";
0087 hPullGM_Y_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0088 title.str("");
0089 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts";
0090 hPullGM_Z_ts[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0091
0092 title.str("");
0093 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr";
0094 hPullGP_X_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0095 title.str("");
0096 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr";
0097 hPullGP_Y_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0098 title.str("");
0099 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr";
0100 hPullGP_Z_tr[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0101
0102 title.str("");
0103 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs";
0104 hPullGP_X_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0105 title.str("");
0106 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs";
0107 hPullGP_Y_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0108 title.str("");
0109 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs";
0110 hPullGP_Z_rs[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0111
0112 if (((i == 2 || i == 4) && (j == 0 || j == 1)) || (i == 3 || i == 5)) {
0113
0114 title.str("");
0115 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
0116 hPullGP_X_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0117 title.str("");
0118 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
0119 hPullGP_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0120 title.str("");
0121 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
0122 hPullGP_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0123
0124 title.str("");
0125 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
0126 hPullGM_X_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0127 title.str("");
0128 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
0129 hPullGM_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0130 title.str("");
0131 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
0132 hPullGM_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0133
0134 title.str("");
0135 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_mono";
0136 hPullGP_X_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0137 title.str("");
0138 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_mono";
0139 hPullGP_Y_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0140 title.str("");
0141 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_mono";
0142 hPullGP_Z_tr_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0143
0144 title.str("");
0145 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_mono";
0146 hPullGP_X_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0147 title.str("");
0148 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_mono";
0149 hPullGP_Y_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0150 title.str("");
0151 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_mono";
0152 hPullGP_Z_rs_mono[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0153
0154
0155 title.str("");
0156 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0157 hPullGP_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0158 title.str("");
0159 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0160 hPullGP_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0161 title.str("");
0162 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0163 hPullGP_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0164
0165 title.str("");
0166 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0167 hPullGM_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0168 title.str("");
0169 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0170 hPullGM_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0171 title.str("");
0172 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0173 hPullGM_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0174
0175 title.str("");
0176 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0177 hPullGP_X_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0178 title.str("");
0179 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0180 hPullGP_Y_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0181 title.str("");
0182 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0183 hPullGP_Z_tr_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0184
0185 title.str("");
0186 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0187 hPullGP_X_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0188 title.str("");
0189 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0190 hPullGP_Y_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0191 title.str("");
0192 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0193 hPullGP_Z_rs_stereo[title.str()] = new TH1F(title.str().c_str(), title.str().c_str(), 1000, -50, 50);
0194 }
0195 }
0196 hTotChi2Increment = new TH1F("TotChi2Increment", "TotChi2Increment", 1000, 0, 100);
0197 hChi2_vs_Process = new TH2F("Chi2_vs_Process", "Chi2_vs_Process", 1000, 0, 100, 17, -0.5, 16.5);
0198 hChi2_vs_clsize = new TH2F("Chi2_vs_clsize", "Chi2_vs_clsize", 1000, 0, 100, 17, -0.5, 16.5);
0199 }
0200
0201 void TestSmoothHits::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0202
0203 edm::ESHandle<TrackerTopology> tTopo = iSetup.getHandle(tTopoToken);
0204
0205 LogTrace("TestSmoothHits") << "new event" << std::endl;
0206
0207 theTCCollection = iEvent.getHandle(theTCCollectionToken);
0208 TrackerHitAssociator hitAssociator(iEvent, trackerHitAssociatorConfig_);
0209
0210 TrajectoryStateCombiner combiner;
0211
0212 for (TrackCandidateCollection::const_iterator i = theTCCollection->begin(); i != theTCCollection->end(); i++) {
0213 LogTrace("TestSmoothHits") << "new candidate" << std::endl;
0214
0215 const TrackCandidate* theTC = &(*i);
0216 PTrajectoryStateOnDet state = theTC->trajectoryStateOnDet();
0217
0218
0219
0220 DetId detId(state.detId());
0221 TrajectoryStateOnSurface theTSOS =
0222 trajectoryStateTransform::transientState(state, &(theG->idToDet(detId)->surface()), theMF.product());
0223
0224 if (theTSOS.globalMomentum().eta() > maxeta || theTSOS.globalMomentum().eta() < mineta)
0225 continue;
0226
0227
0228 TransientTrackingRecHit::RecHitContainer hits;
0229
0230 for (auto const& recHit : theTC->recHits()) {
0231 hits.push_back(theBuilder->build(&recHit));
0232 }
0233
0234
0235 std::vector<Trajectory> fitted = fit->fit(theTC->seed(), hits, theTSOS);
0236
0237 std::vector<Trajectory> result;
0238 for (std::vector<Trajectory>::iterator it = fitted.begin(); it != fitted.end(); it++) {
0239 std::vector<Trajectory> smoothed = smooth->trajectories(*it);
0240 result.insert(result.end(), smoothed.begin(), smoothed.end());
0241 }
0242 if (result.empty())
0243 continue;
0244 std::vector<TrajectoryMeasurement> vtm = result[0].measurements();
0245
0246 TSOS lastState = theTSOS;
0247 for (std::vector<TrajectoryMeasurement>::iterator tm = vtm.begin(); tm != vtm.end(); tm++) {
0248 TransientTrackingRecHit::ConstRecHitPointer rhit = tm->recHit();
0249 if ((rhit)->isValid() == 0 && rhit->det() != nullptr)
0250 continue;
0251 LogTrace("TestSmoothHits") << "new hit";
0252
0253 int subdetId = rhit->det()->geographicalId().subdetId();
0254 DetId id = rhit->det()->geographicalId();
0255 int layerId = tTopo->layer(id);
0256 LogTrace("TestSmoothHits") << "subdetId=" << subdetId << " layerId=" << layerId;
0257
0258 double delta = 99999;
0259 LocalPoint rhitLPv = rhit->localPosition();
0260
0261 std::vector<PSimHit> assSimHits = hitAssociator.associateHit(*(rhit->hit()));
0262 if (assSimHits.empty())
0263 continue;
0264 PSimHit shit;
0265 for (std::vector<PSimHit>::const_iterator m = assSimHits.begin(); m < assSimHits.end(); m++) {
0266 if ((m->localPosition() - rhitLPv).mag() < delta) {
0267 shit = *m;
0268 delta = (m->localPosition() - rhitLPv).mag();
0269 }
0270 }
0271
0272 TSOS currentState = combiner(tm->backwardPredictedState(), tm->forwardPredictedState());
0273 TSOS updatedState = tm->updatedState();
0274
0275
0276 double chi2increment = tm->estimate();
0277 LogTrace("TestSmoothHits") << "tm->estimate()=" << tm->estimate();
0278 title.str("");
0279 title << "Chi2Increment_" << subdetId << "-" << layerId;
0280 hChi2Increment[title.str()]->Fill(chi2increment);
0281 hTotChi2Increment->Fill(chi2increment);
0282 hChi2_vs_Process->Fill(chi2increment, shit.processType());
0283 if (dynamic_cast<const SiPixelRecHit*>(rhit->hit()))
0284 hChi2_vs_clsize->Fill(chi2increment, ((const SiPixelRecHit*)(rhit->hit()))->cluster()->size());
0285 if (dynamic_cast<const SiStripRecHit2D*>(rhit->hit()))
0286 hChi2_vs_clsize->Fill(chi2increment, ((const SiStripRecHit2D*)(rhit->hit()))->cluster()->amplitudes().size());
0287
0288
0289 const Surface* surf = &((rhit)->det()->surface());
0290 LocalVector shitLMom;
0291 LocalPoint shitLPos;
0292 if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())) {
0293 double rechitmatchedx = rhit->localPosition().x();
0294 double rechitmatchedy = rhit->localPosition().y();
0295 double mindist = 999999;
0296 float distx, disty;
0297 std::pair<LocalPoint, LocalVector> closestPair;
0298 const StripGeomDetUnit* stripDet = (StripGeomDetUnit*)((const GluedGeomDet*)(rhit)->det())->stereoDet();
0299 const BoundPlane& plane = (rhit)->det()->surface();
0300 for (std::vector<PSimHit>::const_iterator m = assSimHits.begin(); m < assSimHits.end(); m++) {
0301
0302 std::pair<LocalPoint, LocalVector> hitPair = projectHit((*m), stripDet, plane);
0303 distx = fabs(rechitmatchedx - hitPair.first.x());
0304 disty = fabs(rechitmatchedy - hitPair.first.y());
0305 double dist = distx * distx + disty * disty;
0306 if (sqrt(dist) < mindist) {
0307 mindist = dist;
0308 closestPair = hitPair;
0309 }
0310 }
0311 shitLPos = closestPair.first;
0312 shitLMom = closestPair.second;
0313 } else {
0314 shitLPos = shit.localPosition();
0315 shitLMom = shit.momentumAtEntry();
0316 }
0317 GlobalVector shitGMom = surf->toGlobal(shitLMom);
0318 GlobalPoint shitGPos = surf->toGlobal(shitLPos);
0319
0320
0321
0322
0323
0324
0325
0326
0327
0328
0329
0330
0331
0332
0333
0334
0335
0336
0337
0338
0339
0340
0341
0342
0343
0344
0345
0346
0347
0348
0349
0350
0351
0352
0353
0354
0355
0356
0357
0358
0359
0360 GlobalVector tsosGMom = currentState.globalMomentum();
0361 GlobalError tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
0362 GlobalPoint tsosGPos = currentState.globalPosition();
0363 GlobalError tsosGPEr = currentState.cartesianError().position();
0364
0365 GlobalPoint rhitGPos = (rhit)->globalPosition();
0366 GlobalError rhitGPEr = (rhit)->globalPositionError();
0367
0368 double pullGPX_rs = (rhitGPos.x() - shitGPos.x()) / sqrt(rhitGPEr.cxx());
0369 double pullGPY_rs = (rhitGPos.y() - shitGPos.y()) / sqrt(rhitGPEr.cyy());
0370 double pullGPZ_rs = (rhitGPos.z() - shitGPos.z()) / sqrt(rhitGPEr.czz());
0371
0372
0373
0374
0375 LogTrace("TestSmoothHits") << "rs" << std::endl;
0376
0377 title.str("");
0378 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
0379 hPullGP_X_rs[title.str()]->Fill(pullGPX_rs);
0380 title.str("");
0381 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
0382 hPullGP_Y_rs[title.str()]->Fill(pullGPY_rs);
0383 title.str("");
0384 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
0385 hPullGP_Z_rs[title.str()]->Fill(pullGPZ_rs);
0386
0387 double pullGPX_tr = (tsosGPos.x() - rhitGPos.x()) / sqrt(tsosGPEr.cxx() + rhitGPEr.cxx());
0388 double pullGPY_tr = (tsosGPos.y() - rhitGPos.y()) / sqrt(tsosGPEr.cyy() + rhitGPEr.cyy());
0389 double pullGPZ_tr = (tsosGPos.z() - rhitGPos.z()) / sqrt(tsosGPEr.czz() + rhitGPEr.czz());
0390
0391
0392
0393
0394 LogTrace("TestSmoothHits") << "tr" << std::endl;
0395
0396 title.str("");
0397 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
0398 hPullGP_X_tr[title.str()]->Fill(pullGPX_tr);
0399 title.str("");
0400 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
0401 hPullGP_Y_tr[title.str()]->Fill(pullGPY_tr);
0402 title.str("");
0403 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
0404 hPullGP_Z_tr[title.str()]->Fill(pullGPZ_tr);
0405
0406 double pullGPX_ts = (tsosGPos.x() - shitGPos.x()) / sqrt(tsosGPEr.cxx());
0407 double pullGPY_ts = (tsosGPos.y() - shitGPos.y()) / sqrt(tsosGPEr.cyy());
0408 double pullGPZ_ts = (tsosGPos.z() - shitGPos.z()) / sqrt(tsosGPEr.czz());
0409
0410
0411
0412
0413 LogTrace("TestSmoothHits") << "ts1" << std::endl;
0414
0415 title.str("");
0416 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
0417 hPullGP_X_ts[title.str()]->Fill(pullGPX_ts);
0418 title.str("");
0419 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
0420 hPullGP_Y_ts[title.str()]->Fill(pullGPY_ts);
0421 title.str("");
0422 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
0423 hPullGP_Z_ts[title.str()]->Fill(pullGPZ_ts);
0424
0425 double pullGMX_ts = (tsosGMom.x() - shitGMom.x()) / sqrt(tsosGMEr.cxx());
0426 double pullGMY_ts = (tsosGMom.y() - shitGMom.y()) / sqrt(tsosGMEr.cyy());
0427 double pullGMZ_ts = (tsosGMom.z() - shitGMom.z()) / sqrt(tsosGMEr.czz());
0428
0429
0430
0431
0432 LogTrace("TestSmoothHits") << "ts2" << std::endl;
0433
0434 title.str("");
0435 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
0436 hPullGM_X_ts[title.str()]->Fill(pullGMX_ts);
0437 title.str("");
0438 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
0439 hPullGM_Y_ts[title.str()]->Fill(pullGMY_ts);
0440 title.str("");
0441 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
0442 hPullGM_Z_ts[title.str()]->Fill(pullGMZ_ts);
0443
0444 if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
0445
0446 LogTrace("TestSmoothHits") << "MONO HIT" << std::endl;
0447 auto m = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->monoHit();
0448 CTTRHp tMonoHit = theBuilder->build(&m);
0449 if (tMonoHit == nullptr)
0450 continue;
0451 vector<PSimHit> assMonoSimHits = hitAssociator.associateHit(*tMonoHit->hit());
0452 if (assMonoSimHits.empty())
0453 continue;
0454 const PSimHit sMonoHit = *(assSimHits.begin());
0455 const Surface* monoSurf = &(tMonoHit->det()->surface());
0456 if (monoSurf == nullptr)
0457 continue;
0458 TSOS monoState = thePropagator->propagate(lastState, *monoSurf);
0459 if (monoState.isValid() == 0)
0460 continue;
0461
0462 LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
0463 GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
0464 LocalPoint monoShitLPos = sMonoHit.localPosition();
0465 GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
0466
0467 GlobalVector monoTsosGMom = monoState.globalMomentum();
0468 GlobalError monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
0469 GlobalPoint monoTsosGPos = monoState.globalPosition();
0470 GlobalError monoTsosGPEr = monoState.cartesianError().position();
0471
0472 GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
0473 GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
0474
0475 double pullGPX_rs_mono = (monoRhitGPos.x() - monoShitGPos.x()) / sqrt(monoRhitGPEr.cxx());
0476 double pullGPY_rs_mono = (monoRhitGPos.y() - monoShitGPos.y()) / sqrt(monoRhitGPEr.cyy());
0477 double pullGPZ_rs_mono = (monoRhitGPos.z() - monoShitGPos.z()) / sqrt(monoRhitGPEr.czz());
0478
0479
0480
0481
0482 title.str("");
0483 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
0484 hPullGP_X_rs_mono[title.str()]->Fill(pullGPX_rs_mono);
0485 title.str("");
0486 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
0487 hPullGP_Y_rs_mono[title.str()]->Fill(pullGPY_rs_mono);
0488 title.str("");
0489 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
0490 hPullGP_Z_rs_mono[title.str()]->Fill(pullGPZ_rs_mono);
0491
0492 double pullGPX_tr_mono = (monoTsosGPos.x() - monoRhitGPos.x()) / sqrt(monoTsosGPEr.cxx() + monoRhitGPEr.cxx());
0493 double pullGPY_tr_mono = (monoTsosGPos.y() - monoRhitGPos.y()) / sqrt(monoTsosGPEr.cyy() + monoRhitGPEr.cyy());
0494 double pullGPZ_tr_mono = (monoTsosGPos.z() - monoRhitGPos.z()) / sqrt(monoTsosGPEr.czz() + monoRhitGPEr.czz());
0495
0496
0497
0498
0499 title.str("");
0500 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
0501 hPullGP_X_tr_mono[title.str()]->Fill(pullGPX_tr_mono);
0502 title.str("");
0503 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
0504 hPullGP_Y_tr_mono[title.str()]->Fill(pullGPY_tr_mono);
0505 title.str("");
0506 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
0507 hPullGP_Z_tr_mono[title.str()]->Fill(pullGPZ_tr_mono);
0508
0509 double pullGPX_ts_mono = (monoTsosGPos.x() - monoShitGPos.x()) / sqrt(monoTsosGPEr.cxx());
0510 double pullGPY_ts_mono = (monoTsosGPos.y() - monoShitGPos.y()) / sqrt(monoTsosGPEr.cyy());
0511 double pullGPZ_ts_mono = (monoTsosGPos.z() - monoShitGPos.z()) / sqrt(monoTsosGPEr.czz());
0512
0513
0514
0515
0516 title.str("");
0517 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
0518 hPullGP_X_ts_mono[title.str()]->Fill(pullGPX_ts_mono);
0519 title.str("");
0520 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
0521 hPullGP_Y_ts_mono[title.str()]->Fill(pullGPY_ts_mono);
0522 title.str("");
0523 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
0524 hPullGP_Z_ts_mono[title.str()]->Fill(pullGPZ_ts_mono);
0525
0526 double pullGMX_ts_mono = (monoTsosGMom.x() - monoShitGMom.x()) / sqrt(monoTsosGMEr.cxx());
0527 double pullGMY_ts_mono = (monoTsosGMom.y() - monoShitGMom.y()) / sqrt(monoTsosGMEr.cyy());
0528 double pullGMZ_ts_mono = (monoTsosGMom.z() - monoShitGMom.z()) / sqrt(monoTsosGMEr.czz());
0529
0530
0531
0532
0533 title.str("");
0534 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
0535 hPullGM_X_ts_mono[title.str()]->Fill(pullGMX_ts_mono);
0536 title.str("");
0537 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
0538 hPullGM_Y_ts_mono[title.str()]->Fill(pullGMY_ts_mono);
0539 title.str("");
0540 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
0541 hPullGM_Z_ts_mono[title.str()]->Fill(pullGMZ_ts_mono);
0542
0543
0544 LogTrace("TestSmoothHits") << "STEREO HIT" << std::endl;
0545 auto s = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->stereoHit();
0546 CTTRHp tStereoHit = theBuilder->build(&s);
0547 if (tStereoHit == nullptr)
0548 continue;
0549 vector<PSimHit> assStereoSimHits = hitAssociator.associateHit(*tStereoHit->hit());
0550 if (assStereoSimHits.empty())
0551 continue;
0552 const PSimHit sStereoHit = *(assSimHits.begin());
0553 const Surface* stereoSurf = &(tStereoHit->det()->surface());
0554 if (stereoSurf == nullptr)
0555 continue;
0556 TSOS stereoState = thePropagator->propagate(lastState, *stereoSurf);
0557 if (stereoState.isValid() == 0)
0558 continue;
0559
0560 LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
0561 GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
0562 LocalPoint stereoShitLPos = sStereoHit.localPosition();
0563 GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
0564
0565 GlobalVector stereoTsosGMom = stereoState.globalMomentum();
0566 GlobalError stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3, 3));
0567 GlobalPoint stereoTsosGPos = stereoState.globalPosition();
0568 GlobalError stereoTsosGPEr = stereoState.cartesianError().position();
0569
0570 GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
0571 GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
0572
0573 double pullGPX_rs_stereo = (stereoRhitGPos.x() - stereoShitGPos.x()) / sqrt(stereoRhitGPEr.cxx());
0574 double pullGPY_rs_stereo = (stereoRhitGPos.y() - stereoShitGPos.y()) / sqrt(stereoRhitGPEr.cyy());
0575 double pullGPZ_rs_stereo = (stereoRhitGPos.z() - stereoShitGPos.z()) / sqrt(stereoRhitGPEr.czz());
0576
0577
0578
0579
0580 title.str("");
0581 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
0582 hPullGP_X_rs_stereo[title.str()]->Fill(pullGPX_rs_stereo);
0583 title.str("");
0584 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
0585 hPullGP_Y_rs_stereo[title.str()]->Fill(pullGPY_rs_stereo);
0586 title.str("");
0587 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
0588 hPullGP_Z_rs_stereo[title.str()]->Fill(pullGPZ_rs_stereo);
0589
0590 double pullGPX_tr_stereo =
0591 (stereoTsosGPos.x() - stereoRhitGPos.x()) / sqrt(stereoTsosGPEr.cxx() + stereoRhitGPEr.cxx());
0592 double pullGPY_tr_stereo =
0593 (stereoTsosGPos.y() - stereoRhitGPos.y()) / sqrt(stereoTsosGPEr.cyy() + stereoRhitGPEr.cyy());
0594 double pullGPZ_tr_stereo =
0595 (stereoTsosGPos.z() - stereoRhitGPos.z()) / sqrt(stereoTsosGPEr.czz() + stereoRhitGPEr.czz());
0596
0597
0598
0599
0600 title.str("");
0601 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
0602 hPullGP_X_tr_stereo[title.str()]->Fill(pullGPX_tr_stereo);
0603 title.str("");
0604 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
0605 hPullGP_Y_tr_stereo[title.str()]->Fill(pullGPY_tr_stereo);
0606 title.str("");
0607 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
0608 hPullGP_Z_tr_stereo[title.str()]->Fill(pullGPZ_tr_stereo);
0609
0610 double pullGPX_ts_stereo = (stereoTsosGPos.x() - stereoShitGPos.x()) / sqrt(stereoTsosGPEr.cxx());
0611 double pullGPY_ts_stereo = (stereoTsosGPos.y() - stereoShitGPos.y()) / sqrt(stereoTsosGPEr.cyy());
0612 double pullGPZ_ts_stereo = (stereoTsosGPos.z() - stereoShitGPos.z()) / sqrt(stereoTsosGPEr.czz());
0613
0614
0615
0616
0617 title.str("");
0618 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
0619 hPullGP_X_ts_stereo[title.str()]->Fill(pullGPX_ts_stereo);
0620 title.str("");
0621 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
0622 hPullGP_Y_ts_stereo[title.str()]->Fill(pullGPY_ts_stereo);
0623 title.str("");
0624 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
0625 hPullGP_Z_ts_stereo[title.str()]->Fill(pullGPZ_ts_stereo);
0626
0627 double pullGMX_ts_stereo = (stereoTsosGMom.x() - stereoShitGMom.x()) / sqrt(stereoTsosGMEr.cxx());
0628 double pullGMY_ts_stereo = (stereoTsosGMom.y() - stereoShitGMom.y()) / sqrt(stereoTsosGMEr.cyy());
0629 double pullGMZ_ts_stereo = (stereoTsosGMom.z() - stereoShitGMom.z()) / sqrt(stereoTsosGMEr.czz());
0630
0631
0632
0633
0634 title.str("");
0635 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
0636 hPullGM_X_ts_stereo[title.str()]->Fill(pullGMX_ts_stereo);
0637 title.str("");
0638 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
0639 hPullGM_Y_ts_stereo[title.str()]->Fill(pullGMY_ts_stereo);
0640 title.str("");
0641 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
0642 hPullGM_Z_ts_stereo[title.str()]->Fill(pullGMZ_ts_stereo);
0643 }
0644 lastState = updatedState;
0645
0646 }
0647 }
0648 LogTrace("TestSmoothHits") << "end of event" << std::endl;
0649 }
0650
0651 void TestSmoothHits::endJob() {
0652
0653 TDirectory* chi2i = file->mkdir("Chi2_Increment");
0654
0655 TDirectory* gp_ts = file->mkdir("GP_TSOS-SimHit");
0656 TDirectory* gm_ts = file->mkdir("GM_TSOS-SimHit");
0657 TDirectory* gp_tr = file->mkdir("GP_TSOS-RecHit");
0658 TDirectory* gp_rs = file->mkdir("GP_RecHit-SimHit");
0659
0660 TDirectory* gp_tsx = gp_ts->mkdir("X");
0661 TDirectory* gp_tsy = gp_ts->mkdir("Y");
0662 TDirectory* gp_tsz = gp_ts->mkdir("Z");
0663 TDirectory* gm_tsx = gm_ts->mkdir("X");
0664 TDirectory* gm_tsy = gm_ts->mkdir("Y");
0665 TDirectory* gm_tsz = gm_ts->mkdir("Z");
0666 TDirectory* gp_trx = gp_tr->mkdir("X");
0667 TDirectory* gp_try = gp_tr->mkdir("Y");
0668 TDirectory* gp_trz = gp_tr->mkdir("Z");
0669 TDirectory* gp_rsx = gp_rs->mkdir("X");
0670 TDirectory* gp_rsy = gp_rs->mkdir("Y");
0671 TDirectory* gp_rsz = gp_rs->mkdir("Z");
0672
0673 TDirectory* gp_tsx_mono = gp_ts->mkdir("MONOX");
0674 TDirectory* gp_tsy_mono = gp_ts->mkdir("MONOY");
0675 TDirectory* gp_tsz_mono = gp_ts->mkdir("MONOZ");
0676 TDirectory* gm_tsx_mono = gm_ts->mkdir("MONOX");
0677 TDirectory* gm_tsy_mono = gm_ts->mkdir("MONOY");
0678 TDirectory* gm_tsz_mono = gm_ts->mkdir("MONOZ");
0679 TDirectory* gp_trx_mono = gp_tr->mkdir("MONOX");
0680 TDirectory* gp_try_mono = gp_tr->mkdir("MONOY");
0681 TDirectory* gp_trz_mono = gp_tr->mkdir("MONOZ");
0682 TDirectory* gp_rsx_mono = gp_rs->mkdir("MONOX");
0683 TDirectory* gp_rsy_mono = gp_rs->mkdir("MONOY");
0684 TDirectory* gp_rsz_mono = gp_rs->mkdir("MONOZ");
0685
0686 TDirectory* gp_tsx_stereo = gp_ts->mkdir("STEREOX");
0687 TDirectory* gp_tsy_stereo = gp_ts->mkdir("STEREOY");
0688 TDirectory* gp_tsz_stereo = gp_ts->mkdir("STEREOZ");
0689 TDirectory* gm_tsx_stereo = gm_ts->mkdir("STEREOX");
0690 TDirectory* gm_tsy_stereo = gm_ts->mkdir("STEREOY");
0691 TDirectory* gm_tsz_stereo = gm_ts->mkdir("STEREOZ");
0692 TDirectory* gp_trx_stereo = gp_tr->mkdir("STEREOX");
0693 TDirectory* gp_try_stereo = gp_tr->mkdir("STEREOY");
0694 TDirectory* gp_trz_stereo = gp_tr->mkdir("STEREOZ");
0695 TDirectory* gp_rsx_stereo = gp_rs->mkdir("STEREOX");
0696 TDirectory* gp_rsy_stereo = gp_rs->mkdir("STEREOY");
0697 TDirectory* gp_rsz_stereo = gp_rs->mkdir("STEREOZ");
0698
0699 chi2i->cd();
0700 hTotChi2Increment->Write();
0701 hChi2_vs_Process->Write();
0702 hChi2_vs_clsize->Write();
0703 for (int i = 0; i != 6; i++)
0704 for (int j = 0; j != 9; j++) {
0705 if (i == 0 && j > 2)
0706 break;
0707 if (i == 1 && j > 1)
0708 break;
0709 if (i == 2 && j > 3)
0710 break;
0711 if (i == 3 && j > 2)
0712 break;
0713 if (i == 4 && j > 5)
0714 break;
0715 if (i == 5 && j > 8)
0716 break;
0717 chi2i->cd();
0718 title.str("");
0719 title << "Chi2Increment_" << i + 1 << "-" << j + 1;
0720 hChi2Increment[title.str()]->Write();
0721
0722 gp_ts->cd();
0723 gp_tsx->cd();
0724 title.str("");
0725 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts";
0726 hPullGP_X_ts[title.str()]->Write();
0727 gp_tsy->cd();
0728 title.str("");
0729 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts";
0730 hPullGP_Y_ts[title.str()]->Write();
0731 gp_tsz->cd();
0732 title.str("");
0733 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts";
0734 hPullGP_Z_ts[title.str()]->Write();
0735
0736 gm_ts->cd();
0737 gm_tsx->cd();
0738 title.str("");
0739 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts";
0740 hPullGM_X_ts[title.str()]->Write();
0741 gm_tsy->cd();
0742 title.str("");
0743 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts";
0744 hPullGM_Y_ts[title.str()]->Write();
0745 gm_tsz->cd();
0746 title.str("");
0747 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts";
0748 hPullGM_Z_ts[title.str()]->Write();
0749
0750 gp_tr->cd();
0751 gp_trx->cd();
0752 title.str("");
0753 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr";
0754 hPullGP_X_tr[title.str()]->Write();
0755 gp_try->cd();
0756 title.str("");
0757 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr";
0758 hPullGP_Y_tr[title.str()]->Write();
0759 gp_trz->cd();
0760 title.str("");
0761 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr";
0762 hPullGP_Z_tr[title.str()]->Write();
0763
0764 gp_rs->cd();
0765 gp_rsx->cd();
0766 title.str("");
0767 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs";
0768 hPullGP_X_rs[title.str()]->Write();
0769 gp_rsy->cd();
0770 title.str("");
0771 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs";
0772 hPullGP_Y_rs[title.str()]->Write();
0773 gp_rsz->cd();
0774 title.str("");
0775 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs";
0776 hPullGP_Z_rs[title.str()]->Write();
0777
0778 if (((i == 2 || i == 4) && (j == 0 || j == 1)) || (i == 3 || i == 5)) {
0779
0780 gp_ts->cd();
0781 gp_tsx_mono->cd();
0782 title.str("");
0783 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
0784 hPullGP_X_ts_mono[title.str()]->Write();
0785 gp_tsy_mono->cd();
0786 title.str("");
0787 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
0788 hPullGP_Y_ts_mono[title.str()]->Write();
0789 gp_tsz_mono->cd();
0790 title.str("");
0791 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
0792 hPullGP_Z_ts_mono[title.str()]->Write();
0793
0794 gm_ts->cd();
0795 gm_tsx_mono->cd();
0796 title.str("");
0797 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_mono";
0798 hPullGM_X_ts_mono[title.str()]->Write();
0799 gm_tsy_mono->cd();
0800 title.str("");
0801 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_mono";
0802 hPullGM_Y_ts_mono[title.str()]->Write();
0803 gm_tsz_mono->cd();
0804 title.str("");
0805 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_mono";
0806 hPullGM_Z_ts_mono[title.str()]->Write();
0807
0808 gp_tr->cd();
0809 gp_trx_mono->cd();
0810 title.str("");
0811 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_mono";
0812 hPullGP_X_tr_mono[title.str()]->Write();
0813 gp_try_mono->cd();
0814 title.str("");
0815 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_mono";
0816 hPullGP_Y_tr_mono[title.str()]->Write();
0817 gp_trz_mono->cd();
0818 title.str("");
0819 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_mono";
0820 hPullGP_Z_tr_mono[title.str()]->Write();
0821
0822 gp_rs->cd();
0823 gp_rsx_mono->cd();
0824 title.str("");
0825 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_mono";
0826 hPullGP_X_rs_mono[title.str()]->Write();
0827 gp_rsy_mono->cd();
0828 title.str("");
0829 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_mono";
0830 hPullGP_Y_rs_mono[title.str()]->Write();
0831 gp_rsz_mono->cd();
0832 title.str("");
0833 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_mono";
0834 hPullGP_Z_rs_mono[title.str()]->Write();
0835
0836
0837 gp_ts->cd();
0838 gp_tsx_stereo->cd();
0839 title.str("");
0840 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0841 hPullGP_X_ts_stereo[title.str()]->Write();
0842 gp_tsy_stereo->cd();
0843 title.str("");
0844 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0845 hPullGP_Y_ts_stereo[title.str()]->Write();
0846 gp_tsz_stereo->cd();
0847 title.str("");
0848 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0849 hPullGP_Z_ts_stereo[title.str()]->Write();
0850
0851 gm_ts->cd();
0852 gm_tsx_stereo->cd();
0853 title.str("");
0854 title << "PullGM_X_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0855 hPullGM_X_ts_stereo[title.str()]->Write();
0856 gm_tsy_stereo->cd();
0857 title.str("");
0858 title << "PullGM_Y_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0859 hPullGM_Y_ts_stereo[title.str()]->Write();
0860 gm_tsz_stereo->cd();
0861 title.str("");
0862 title << "PullGM_Z_" << i + 1 << "-" << j + 1 << "_ts_stereo";
0863 hPullGM_Z_ts_stereo[title.str()]->Write();
0864
0865 gp_tr->cd();
0866 gp_trx_stereo->cd();
0867 title.str("");
0868 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0869 hPullGP_X_tr_stereo[title.str()]->Write();
0870 gp_try_stereo->cd();
0871 title.str("");
0872 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0873 hPullGP_Y_tr_stereo[title.str()]->Write();
0874 gp_trz_stereo->cd();
0875 title.str("");
0876 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_tr_stereo";
0877 hPullGP_Z_tr_stereo[title.str()]->Write();
0878
0879 gp_rs->cd();
0880 gp_rsx_stereo->cd();
0881 title.str("");
0882 title << "PullGP_X_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0883 hPullGP_X_rs_stereo[title.str()]->Write();
0884 gp_rsy_stereo->cd();
0885 title.str("");
0886 title << "PullGP_Y_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0887 hPullGP_Y_rs_stereo[title.str()]->Write();
0888 gp_rsz_stereo->cd();
0889 title.str("");
0890 title << "PullGP_Z_" << i + 1 << "-" << j + 1 << "_rs_stereo";
0891 hPullGP_Z_rs_stereo[title.str()]->Write();
0892 }
0893 }
0894
0895 file->Close();
0896 }
0897
0898
0899
0900 std::pair<LocalPoint, LocalVector> TestSmoothHits::projectHit(const PSimHit& hit,
0901 const StripGeomDetUnit* stripDet,
0902 const BoundPlane& plane) {
0903 const StripTopology& topol = stripDet->specificTopology();
0904 GlobalPoint globalpos = stripDet->surface().toGlobal(hit.localPosition());
0905 LocalPoint localHit = plane.toLocal(globalpos);
0906
0907 LocalVector locdir = hit.localDirection();
0908
0909
0910 GlobalVector globaldir = stripDet->surface().toGlobal(locdir);
0911 LocalVector dir = plane.toLocal(globaldir);
0912 float scale = -localHit.z() / dir.z();
0913
0914 LocalPoint projectedPos = localHit + scale * dir;
0915
0916 float selfAngle = topol.stripAngle(topol.strip(hit.localPosition()));
0917
0918 LocalVector stripDir(sin(selfAngle), cos(selfAngle), 0);
0919
0920 LocalVector localStripDir(plane.toLocal(stripDet->surface().toGlobal(stripDir)));
0921
0922 return std::pair<LocalPoint, LocalVector>(projectedPos, localStripDir);
0923 }
0924
0925 #include "FWCore/Framework/interface/ModuleFactory.h"
0926 #include "FWCore/Framework/interface/MakerMacros.h"
0927
0928 DEFINE_FWK_MODULE(TestSmoothHits);
