Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2023-03-17 11:22:11 

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

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!