Reference Manual

Go to the documentation of this file.
00001 #include "RecoTracker/DebugTools/interface/TestHits.h"
00002 
00003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00004 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00005 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
00006 #include "TrackingTools/PatternTools/interface/TransverseImpactPointExtrapolator.h"
00007 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00008 #include "RecoTracker/TrackProducer/interface/TrackingRecHitLessFromGlobalPosition.h"
00009 #include "TrackingTools/PatternTools/interface/TSCPBuilderNoMaterial.h"
00010 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
00011 #include <TDirectory.h>
00012 #include "Geometry/TrackerGeometryBuilder/interface/GluedGeomDet.h"
00013 #include "TrackingTools/TrackFitters/interface/KFTrajectoryFitter.h"
00014 
00015 #include "TrackingTools/TrackFitters/interface/TrajectoryFitter.h"
00016 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
00017 #include "Geometry/Records/interface/IdealGeometryRecord.h"
00018 
00019 typedef TrajectoryStateOnSurface TSOS;
00020 typedef TransientTrackingRecHit::ConstRecHitPointer CTTRHp;
00021 using namespace std;
00022 using namespace edm;
00023 
00024 TestHits::TestHits(const edm::ParameterSet& iConfig):
00025   conf_(iConfig){
00026   LogTrace("TestHits") << conf_<< std::endl;
00027   propagatorName = conf_.getParameter<std::string>("Propagator");   
00028   builderName = conf_.getParameter<std::string>("TTRHBuilder");   
00029   srcName = conf_.getParameter<std::string>("src");   
00030   fname = conf_.getParameter<std::string>("Fitter");
00031   mineta = conf_.getParameter<double>("mineta");
00032   maxeta = conf_.getParameter<double>("maxeta");
00033 }
00034 
00035 TestHits::~TestHits(){}
00036 
00037 void TestHits::beginRun(edm::Run & run, const edm::EventSetup& iSetup)
00038 {
00039  
00040   iSetup.get<TrackerDigiGeometryRecord>().get(theG);
00041   iSetup.get<IdealMagneticFieldRecord>().get(theMF);  
00042   iSetup.get<TrackingComponentsRecord>().get(propagatorName,thePropagator);
00043   iSetup.get<TransientRecHitRecord>().get(builderName,theBuilder);
00044   iSetup.get<TrajectoryFitter::Record>().get(fname, fit);
00045  
00046   file = new TFile("testhits.root","recreate");
00047   for (int i=0; i!=6; i++)
00048     for (int j=0; j!=9; j++){
00049       if (i==0 && j>2) break;
00050       if (i==1 && j>1) break;
00051       if (i==2 && j>3) break;
00052       if (i==3 && j>2) break;
00053       if (i==4 && j>5) break;
00054       if (i==5 && j>8) break;
00055       title.str("");
00056       title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
00057       hPullGP_X_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00058       title.str("");
00059       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts";
00060       hPullGP_Y_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00061       title.str("");
00062       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts";
00063       hPullGP_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00064       title.str("");
00065       title << "Chi2Increment_" << i+1 << "-" << j+1;
00066       hChi2Increment[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00067 
00068       title.str("");
00069       title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
00070       hPullGM_X_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00071       title.str("");
00072       title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
00073       hPullGM_Y_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00074       title.str("");
00075       title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
00076       hPullGM_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00077 
00078       title.str("");
00079       title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
00080       hPullGP_X_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00081       title.str("");
00082       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
00083       hPullGP_Y_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00084       title.str("");
00085       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
00086       hPullGP_Z_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00087 
00088       title.str("");
00089       title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
00090       hPullGP_X_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00091       title.str("");
00092       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
00093       hPullGP_Y_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00094       title.str("");
00095       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
00096       hPullGP_Z_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00097 
00098       if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
00099         //mono
00100         title.str("");
00101         title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
00102         hPullGP_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00103         title.str("");
00104         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
00105         hPullGP_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00106         title.str("");
00107         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
00108         hPullGP_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00109 
00110         title.str("");
00111         title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
00112         hPullGM_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00113         title.str("");
00114         title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
00115         hPullGM_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00116         title.str("");
00117         title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
00118         hPullGM_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00119 
00120         title.str("");
00121         title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
00122         hPullGP_X_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00123         title.str("");
00124         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
00125         hPullGP_Y_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00126         title.str("");
00127         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
00128         hPullGP_Z_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00129 
00130         title.str("");
00131         title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
00132         hPullGP_X_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00133         title.str("");
00134         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
00135         hPullGP_Y_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00136         title.str("");
00137         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
00138         hPullGP_Z_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00139 
00140         //stereo
00141         title.str("");
00142         title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
00143         hPullGP_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00144         title.str("");
00145         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
00146         hPullGP_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00147         title.str("");
00148         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
00149         hPullGP_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00150 
00151         title.str("");
00152         title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
00153         hPullGM_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00154         title.str("");
00155         title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
00156         hPullGM_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00157         title.str("");
00158         title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
00159         hPullGM_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00160 
00161         title.str("");
00162         title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
00163         hPullGP_X_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00164         title.str("");
00165         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
00166         hPullGP_Y_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00167         title.str("");
00168         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
00169         hPullGP_Z_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00170 
00171         title.str("");
00172         title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
00173         hPullGP_X_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00174         title.str("");
00175         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
00176         hPullGP_Y_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00177         title.str("");
00178         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
00179         hPullGP_Z_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00180       }
00181     }
00182   hTotChi2Increment = new TH1F("TotChi2Increment","TotChi2Increment",1000,0,100);
00183   hProcess_vs_Chi2  = new TH2F("Process_vs_Chi2","Process_vs_Chi2",1000,0,100,17,-0.5,16.5);  
00184   hClsize_vs_Chi2  = new TH2F("Clsize_vs_Chi2","Clsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00185 }
00186 
00187 
00188 void TestHits::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00189 {
00190   //Retrieve tracker topology from geometry
00191   edm::ESHandle<TrackerTopology> tTopo;
00192   iSetup.get<IdealGeometryRecord>().get(tTopo);
00193 
00194 
00195   LogTrace("TestHits") << "\nnew event";
00196 
00197   iEvent.getByLabel(srcName,theTCCollection ); 
00198   hitAssociator = new TrackerHitAssociator(iEvent);
00199 
00200   TrajectoryStateCombiner combiner;
00201 
00202   for (TrackCandidateCollection::const_iterator i=theTCCollection->begin(); i!=theTCCollection->end();i++){
00203 
00204     LogTrace("TestHits") << "\n*****************new candidate*****************" << std::endl;
00205       
00206     const TrackCandidate * theTC = &(*i);
00207     PTrajectoryStateOnDet state = theTC->trajectoryStateOnDet();
00208     const TrackCandidate::range& recHitVec=theTC->recHits();
00209 
00210     //convert PTrajectoryStateOnDet to TrajectoryStateOnSurface
00211     
00212     
00213     DetId  detId(state.detId());
00214     TrajectoryStateOnSurface theTSOS=
00215       trajectoryStateTransform::transientState(state, &(theG->idToDet(detId)->surface()),theMF.product());
00216 
00217     if (theTSOS.globalMomentum().eta()>maxeta || theTSOS.globalMomentum().eta()<mineta) continue;
00218     
00219     //convert the TrackingRecHit vector to a TransientTrackingRecHit vector
00220     TransientTrackingRecHit::RecHitContainer hits;
00221     
00222     for (edm::OwnVector<TrackingRecHit>::const_iterator i=recHitVec.first;
00223          i!=recHitVec.second; i++){
00224       hits.push_back(theBuilder->build(&(*i) ));
00225     }
00226 
00227     //call the fitter
00228     std::vector<Trajectory> result = fit->fit(theTC->seed(), hits, theTSOS);
00229     if (result.size()==0) continue;
00230     std::vector<TrajectoryMeasurement> vtm = result[0].measurements();
00231     double tchi2 = 0;
00232 
00233     TSOS lastState = theTSOS;
00234     for (std::vector<TrajectoryMeasurement>::iterator tm=vtm.begin(); tm!=vtm.end();tm++){
00235 
00236       TransientTrackingRecHit::ConstRecHitPointer rhit = tm->recHit();
00237       if ((rhit)->isValid()==0&&rhit->det()!=0) continue;
00238       LogTrace("TestHits") << "*****************new hit*****************" ;
00239 
00240       int subdetId = rhit->det()->geographicalId().subdetId();
00241       DetId id = rhit->det()->geographicalId();
00242       int layerId  = tTopo->layer(id);
00243       LogTrace("TestHits") << "subdetId=" << subdetId << " layerId=" << layerId ;
00244 
00245       double delta = 99999;
00246       LocalPoint rhitLPv = rhit->localPosition();
00247 
00248       std::vector<PSimHit> assSimHits = hitAssociator->associateHit(*(rhit->hit()));
00249       if (assSimHits.size()==0) continue;
00250       PSimHit shit;
00251       for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
00252         if ((m->localPosition()-rhitLPv).mag()<delta) {
00253           shit=*m;
00254           delta = (m->localPosition()-rhitLPv).mag();
00255         }
00256       }
00257 
00258       TSOS currentState = tm->forwardPredictedState();
00259       if (tm->backwardPredictedState().isValid()) 
00260         currentState = combiner(tm->backwardPredictedState(), tm->forwardPredictedState());
00261       TSOS updatedState = tm->updatedState();
00262       tchi2+=tm->estimate();
00263 
00264       //plot chi2 increment
00265       double chi2increment = tm->estimate();
00266       LogTrace("TestHits") << "tm->estimate()=" << tm->estimate();
00267       title.str("");
00268       title << "Chi2Increment_" << subdetId << "-" << layerId;
00269       hChi2Increment[title.str()]->Fill( chi2increment );
00270       hTotChi2Increment->Fill( chi2increment );
00271       hProcess_vs_Chi2->Fill( chi2increment, shit.processType() );
00272       if (dynamic_cast<const SiPixelRecHit*>(rhit->hit()))      
00273         hClsize_vs_Chi2->Fill( chi2increment, ((const SiPixelRecHit*)(rhit->hit()))->cluster()->size() );
00274       if (dynamic_cast<const SiStripRecHit2D*>(rhit->hit()))    
00275         hClsize_vs_Chi2->Fill( chi2increment, ((const SiStripRecHit2D*)(rhit->hit()))->cluster()->amplitudes().size() );
00276      
00277       //test hits
00278       const Surface * surf = &( (rhit)->det()->surface() );
00279       LocalVector shitLMom;
00280       LocalPoint shitLPos;
00281       if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
00282         double rechitmatchedx = rhit->localPosition().x();
00283         double rechitmatchedy = rhit->localPosition().y();
00284         double mindist = 999999;
00285         double distx, disty;
00286         std::pair<LocalPoint,LocalVector> closestPair;
00287         const StripGeomDetUnit* stripDet =(StripGeomDetUnit*) ((const GluedGeomDet *)(rhit)->det())->stereoDet();
00288         const BoundPlane& plane = (rhit)->det()->surface();
00289         for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++) {
00290           //project simhit;
00291           std::pair<LocalPoint,LocalVector> hitPair = projectHit((*m),stripDet,plane);
00292           distx = fabs(rechitmatchedx - hitPair.first.x());
00293           disty = fabs(rechitmatchedy - hitPair.first.y());
00294           double dist = distx*distx+disty*disty;
00295           if(sqrt(dist)<mindist){
00296             mindist = dist;
00297             closestPair = hitPair;
00298           }
00299         }
00300         shitLPos = closestPair.first;
00301         shitLMom = closestPair.second;
00302       } else {
00303         shitLPos = shit.localPosition();        
00304         shitLMom = shit.momentumAtEntry();
00305       }
00306       GlobalVector shitGMom = surf->toGlobal(shitLMom);
00307       GlobalPoint shitGPos = surf->toGlobal(shitLPos);
00308 
00309       GlobalVector tsosGMom = currentState.globalMomentum();
00310       GlobalError  tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00311       GlobalPoint  tsosGPos = currentState.globalPosition();
00312       GlobalError  tsosGPEr = currentState.cartesianError().position();
00313 
00314       GlobalPoint rhitGPos = (rhit)->globalPosition();
00315       GlobalError rhitGPEr = (rhit)->globalPositionError();
00316 
00317       double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/sqrt(rhitGPEr.cxx());
00318       double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/sqrt(rhitGPEr.cyy());
00319       double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/sqrt(rhitGPEr.czz());
00320       //double pullGPX_rs = (rhitGPos.x()-shitGPos.x());
00321       //double pullGPY_rs = (rhitGPos.y()-shitGPos.y());
00322       //double pullGPZ_rs = (rhitGPos.z()-shitGPos.z());
00323       
00324       LogTrace("TestHits") << "rs" << std::endl;
00325       LogVerbatim("TestHits") << "assSimHits.size()=" << assSimHits.size() ;
00326       LogVerbatim("TestHits") << "tsos globalPos   =" << tsosGPos ;
00327       LogVerbatim("TestHits") << "sim hit globalPos=" << shitGPos ;
00328       LogVerbatim("TestHits") << "rec hit globalPos=" << rhitGPos ;
00329       LogVerbatim("TestHits") << "geographicalId   =" << rhit->det()->geographicalId().rawId() ;
00330       LogVerbatim("TestHits") << "surface position =" << surf->position() ;
00331 
00332       title.str("");
00333       title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
00334       hPullGP_X_rs[title.str()]->Fill( pullGPX_rs );
00335       title.str("");
00336       title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
00337       hPullGP_Y_rs[title.str()]->Fill( pullGPY_rs );
00338       title.str("");
00339       title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
00340       hPullGP_Z_rs[title.str()]->Fill( pullGPZ_rs );
00341 
00342       double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx());
00343       double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy());
00344       double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/sqrt(tsosGPEr.czz()+rhitGPEr.czz());
00345       //double pullGPX_tr = (tsosGPos.x()-rhitGPos.x());
00346       //double pullGPY_tr = (tsosGPos.y()-rhitGPos.y());
00347       //double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z());
00348 
00349       LogTrace("TestHits") << "tr" << std::endl;
00350 
00351       title.str("");
00352       title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
00353       hPullGP_X_tr[title.str()]->Fill( pullGPX_tr );
00354       title.str("");
00355       title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
00356       hPullGP_Y_tr[title.str()]->Fill( pullGPY_tr );
00357       title.str("");
00358       title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
00359       hPullGP_Z_tr[title.str()]->Fill( pullGPZ_tr );
00360 
00361       double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/sqrt(tsosGPEr.cxx());
00362       double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/sqrt(tsosGPEr.cyy());
00363       double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/sqrt(tsosGPEr.czz());
00364       //double pullGPX_ts = (tsosGPos.x()-shitGPos.x());
00365       //double pullGPY_ts = (tsosGPos.y()-shitGPos.y());
00366       //double pullGPZ_ts = (tsosGPos.z()-shitGPos.z());
00367 
00368       LogTrace("TestHits") << "ts1" << std::endl;
00369 
00370       title.str("");
00371       title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
00372       hPullGP_X_ts[title.str()]->Fill( pullGPX_ts );
00373       title.str("");
00374       title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
00375       hPullGP_Y_ts[title.str()]->Fill( pullGPY_ts );
00376       title.str("");
00377       title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
00378       hPullGP_Z_ts[title.str()]->Fill( pullGPZ_ts );
00379 
00380       double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/sqrt(tsosGMEr.cxx());
00381       double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/sqrt(tsosGMEr.cyy());
00382       double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/sqrt(tsosGMEr.czz());
00383       //double pullGMX_ts = (tsosGMom.x()-shitGMom.x());
00384       //double pullGMY_ts = (tsosGMom.y()-shitGMom.y());
00385       //double pullGMZ_ts = (tsosGMom.z()-shitGMom.z());
00386 
00387       LogTrace("TestHits") << "ts2" << std::endl;
00388 
00389       title.str("");
00390       title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
00391       hPullGM_X_ts[title.str()]->Fill( pullGMX_ts );
00392       title.str("");
00393       title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
00394       hPullGM_Y_ts[title.str()]->Fill( pullGMY_ts );
00395       title.str("");
00396       title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
00397       hPullGM_Z_ts[title.str()]->Fill( pullGMZ_ts );
00398 
00399       if (dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())) {
00400         //mono
00401         LogTrace("TestHits") << "MONO HIT" << std::endl;
00402         auto m = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->monoHit();
00403         CTTRHp tMonoHit = theBuilder->build(&m);
00404         if (tMonoHit==0) continue;
00405         vector<PSimHit> assMonoSimHits = hitAssociator->associateHit(*tMonoHit->hit());
00406         if (assMonoSimHits.size()==0) continue;
00407         const PSimHit sMonoHit = *(assSimHits.begin());
00408         const Surface * monoSurf = &( tMonoHit->det()->surface() );
00409         if (monoSurf==0) continue;
00410         TSOS monoState = thePropagator->propagate(lastState,*monoSurf);
00411         if (monoState.isValid()==0) continue;
00412 
00413         LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
00414         GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
00415         LocalPoint monoShitLPos = sMonoHit.localPosition();
00416         GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
00417 
00418         GlobalVector monoTsosGMom = monoState.globalMomentum();
00419         GlobalError  monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00420         GlobalPoint  monoTsosGPos = monoState.globalPosition();
00421         GlobalError  monoTsosGPEr = monoState.cartesianError().position();
00422 
00423         GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
00424         GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
00425 
00426         double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/sqrt(monoRhitGPEr.cxx());
00427         double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/sqrt(monoRhitGPEr.cyy());
00428         double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/sqrt(monoRhitGPEr.czz());
00429         //double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x());
00430         //double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y());
00431         //double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z());
00432 
00433         title.str("");
00434         title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
00435         hPullGP_X_rs_mono[title.str()]->Fill( pullGPX_rs_mono );
00436         title.str("");
00437         title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
00438         hPullGP_Y_rs_mono[title.str()]->Fill( pullGPY_rs_mono );
00439         title.str("");
00440         title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
00441         hPullGP_Z_rs_mono[title.str()]->Fill( pullGPZ_rs_mono );
00442 
00443         double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx());
00444         double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy());
00445         double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz());
00446         //double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x());
00447         //double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y());
00448         //double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z());
00449 
00450         title.str("");
00451         title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
00452         hPullGP_X_tr_mono[title.str()]->Fill( pullGPX_tr_mono );
00453         title.str("");
00454         title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
00455         hPullGP_Y_tr_mono[title.str()]->Fill( pullGPY_tr_mono );
00456         title.str("");
00457         title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
00458         hPullGP_Z_tr_mono[title.str()]->Fill( pullGPZ_tr_mono );
00459 
00460         double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/sqrt(monoTsosGPEr.cxx());
00461         double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/sqrt(monoTsosGPEr.cyy());
00462         double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/sqrt(monoTsosGPEr.czz());
00463         //double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x());
00464         //double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y());
00465         //double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z());
00466 
00467         title.str("");
00468         title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
00469         hPullGP_X_ts_mono[title.str()]->Fill( pullGPX_ts_mono );
00470         title.str("");
00471         title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
00472         hPullGP_Y_ts_mono[title.str()]->Fill( pullGPY_ts_mono );
00473         title.str("");
00474         title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
00475         hPullGP_Z_ts_mono[title.str()]->Fill( pullGPZ_ts_mono );
00476 
00477         double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/sqrt(monoTsosGMEr.cxx());
00478         double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/sqrt(monoTsosGMEr.cyy());
00479         double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/sqrt(monoTsosGMEr.czz());
00480         //double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x());
00481         //double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y());
00482         //double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z());
00483 
00484         title.str("");
00485         title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
00486         hPullGM_X_ts_mono[title.str()]->Fill( pullGMX_ts_mono );
00487         title.str("");
00488         title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
00489         hPullGM_Y_ts_mono[title.str()]->Fill( pullGMY_ts_mono );
00490         title.str("");
00491         title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
00492         hPullGM_Z_ts_mono[title.str()]->Fill( pullGMZ_ts_mono );
00493 
00494         //stereo
00495         LogTrace("TestHits") << "STEREO HIT" << std::endl;
00496         auto s = dynamic_cast<const SiStripMatchedRecHit2D*>((rhit)->hit())->stereoHit();
00497         CTTRHp tStereoHit = 
00498           theBuilder->build(&s);
00499         if (tStereoHit==0) continue;
00500         vector<PSimHit> assStereoSimHits = hitAssociator->associateHit(*tStereoHit->hit());
00501         if (assStereoSimHits.size()==0) continue;
00502         const PSimHit sStereoHit = *(assSimHits.begin());
00503         const Surface * stereoSurf = &( tStereoHit->det()->surface() );
00504         if (stereoSurf==0) continue;
00505         TSOS stereoState = thePropagator->propagate(lastState,*stereoSurf);
00506         if (stereoState.isValid()==0) continue;
00507 
00508         LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
00509         GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
00510         LocalPoint stereoShitLPos = sStereoHit.localPosition();
00511         GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
00512 
00513         GlobalVector stereoTsosGMom = stereoState.globalMomentum();
00514         GlobalError  stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00515         GlobalPoint  stereoTsosGPos = stereoState.globalPosition();
00516         GlobalError  stereoTsosGPEr = stereoState.cartesianError().position();
00517 
00518         GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
00519         GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
00520 
00521         double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/sqrt(stereoRhitGPEr.cxx());
00522         double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/sqrt(stereoRhitGPEr.cyy());
00523         double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/sqrt(stereoRhitGPEr.czz());
00524         //double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x());
00525         //double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y());
00526         //double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z());
00527 
00528         title.str("");
00529         title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
00530         hPullGP_X_rs_stereo[title.str()]->Fill( pullGPX_rs_stereo );
00531         title.str("");
00532         title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
00533         hPullGP_Y_rs_stereo[title.str()]->Fill( pullGPY_rs_stereo );
00534         title.str("");
00535         title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
00536         hPullGP_Z_rs_stereo[title.str()]->Fill( pullGPZ_rs_stereo );
00537 
00538         double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx());
00539         double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy());
00540         double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz());
00541         //double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x());
00542         //double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y());
00543         //double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z());
00544 
00545         title.str("");
00546         title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
00547         hPullGP_X_tr_stereo[title.str()]->Fill( pullGPX_tr_stereo );
00548         title.str("");
00549         title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
00550         hPullGP_Y_tr_stereo[title.str()]->Fill( pullGPY_tr_stereo );
00551         title.str("");
00552         title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
00553         hPullGP_Z_tr_stereo[title.str()]->Fill( pullGPZ_tr_stereo );
00554 
00555         double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/sqrt(stereoTsosGPEr.cxx());
00556         double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/sqrt(stereoTsosGPEr.cyy());
00557         double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/sqrt(stereoTsosGPEr.czz());
00558         //double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x());
00559         //double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y());
00560         //double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z());
00561 
00562         title.str("");
00563         title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
00564         hPullGP_X_ts_stereo[title.str()]->Fill( pullGPX_ts_stereo );
00565         title.str("");
00566         title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00567         hPullGP_Y_ts_stereo[title.str()]->Fill( pullGPY_ts_stereo );
00568         title.str("");
00569         title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00570         hPullGP_Z_ts_stereo[title.str()]->Fill( pullGPZ_ts_stereo );
00571 
00572         double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/sqrt(stereoTsosGMEr.cxx());
00573         double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/sqrt(stereoTsosGMEr.cyy());
00574         double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/sqrt(stereoTsosGMEr.czz());
00575         //double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x());
00576         //double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y());
00577         //double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z());
00578 
00579         title.str("");
00580         title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
00581         hPullGM_X_ts_stereo[title.str()]->Fill( pullGMX_ts_stereo );
00582         title.str("");
00583         title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00584         hPullGM_Y_ts_stereo[title.str()]->Fill( pullGMY_ts_stereo );
00585         title.str("");
00586         title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00587         hPullGM_Z_ts_stereo[title.str()]->Fill( pullGMZ_ts_stereo );
00588       }    
00589       lastState = updatedState;
00590     }
00591     LogTrace("TestHits") << "traj chi2="  << tchi2 ;
00592     LogTrace("TestHits") << "track chi2=" << result[0].chiSquared() ;
00593   }
00594   delete hitAssociator;
00595   LogTrace("TestHits") << "end of event" << std::endl;
00596 }
00597 //     TSOS lastState = theTSOS;
00598 //     for (std::vector<TrajectoryMeasurement>::iterator tm=vtm.begin(); tm!=vtm.end();tm++){
00599 
00600 //       TransientTrackingRecHit::ConstRecHitPointer rhit = tm->recHit();
00601 //       if ((rhit)->isValid()==0&&rhit->det()!=0) continue;
00602     
00603 // //     //test hits
00604 // //     TSOS lastState, currentState, updatedState;
00605 // //     updatedState=theTSOS;
00606 // //     for (TransientTrackingRecHit::RecHitContainer::iterator rhit=hits.begin()+1;rhit!=hits.end();++rhit){
00607 
00608 //       lastState=updatedState;
00609 
00610 //       LogTrace("TestHits") << "new hit" << std::endl;
00611 
00612 //       if ((*rhit)->isValid()==0) continue;
00613 
00614 //       int subdetId = (*rhit)->det()->geographicalId().subdetId();
00615 //       int layerId  = 0;
00616 //       DetId id = (*rhit)->det()->geographicalId();
00617 //       if (id.subdetId()==3) layerId = ((tTopo->tibLayer(id);
00618 //       if (id.subdetId()==5) layerId = ((tTopo->tobLayer(id);
00619 //       if (id.subdetId()==1) layerId = ((tTopo->pxbLayer(id);
00620 //       if (id.subdetId()==4) layerId = ((tTopo->tidWheel(id);
00621 //       if (id.subdetId()==6) layerId = ((tTopo->tecWheel(id);
00622 //       if (id.subdetId()==2) layerId = ((tTopo->pxfDisk(id);
00623 //       const Surface * surf = &( (*rhit)->det()->surface() );
00624 //       currentState=thePropagator->propagate(lastState,*surf);        
00625 //       if (currentState.isValid()==0) continue;
00626 //       updatedState=theUpdator->update(currentState,**rhit);
00627 
00628 //       double delta = 99999;
00629 //       LocalPoint rhitLP = rhit->localPosition();
00630 
00631 //       std::vector<PSimHit> assSimHits = hitAssociator->associateHit(*(*rhit)->hit());
00632 //       if (assSimHits.size()==0) continue;
00633 //       PSimHit shit;
00634 //       for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
00635 //      if ((*m-rhitLP).mag()<delta) {
00636 //        shit=*m;
00637 //        delta = (*m-rhitLP).mag();
00638 //      }
00639 //       }
00640 
00641 //       LocalVector shitLMom = shit.momentumAtEntry();
00642 //       GlobalVector shitGMom = surf->toGlobal(shitLMom);
00643 //       LocalPoint shitLPos = shit.localPosition();
00644 //       GlobalPoint shitGPos = surf->toGlobal(shitLPos);
00645 
00646 //       GlobalVector tsosGMom = currentState.globalMomentum();
00647 //       GlobalError  tsosGMEr(currentState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00648 //       GlobalPoint  tsosGPos = currentState.globalPosition();
00649 //       GlobalError  tsosGPEr = currentState.cartesianError().position();
00650 
00651 //       GlobalPoint rhitGPos = (*rhit)->globalPosition();
00652 //       GlobalError rhitGPEr = (*rhit)->globalPositionError();
00653 
00654 //       double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/sqrt(rhitGPEr.cxx());
00655 //       double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/sqrt(rhitGPEr.cyy());
00656 //       double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/sqrt(rhitGPEr.czz());
00657 // //       double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/*/sqrt(rhitGPEr.cxx())*/;
00658 // //       double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/*/sqrt(rhitGPEr.cyy())*/;
00659 // //       double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/*/sqrt(rhitGPEr.czz())*/;
00660 
00661 //       //plot chi2 increment
00662 //       MeasurementExtractor me(currentState);
00663 //       double chi2increment = computeChi2Increment(me,*rhit);
00664 //       LogTrace("TestHits") << "chi2increment=" << chi2increment << std::endl;
00665 //       title.str("");
00666 //       title << "Chi2Increment_" << subdetId << "-" << layerId;
00667 //       hChi2Increment[title.str()]->Fill( chi2increment );
00668 //       hTotChi2Increment->Fill( chi2increment );
00669 //       hChi2_vs_Process->Fill( chi2increment, shit.processType() );
00670 //       if (dynamic_cast<const SiPixelRecHit*>((*rhit)->hit()))        
00671 //      hChi2_vs_clsize->Fill( chi2increment, ((const SiPixelRecHit*)(*rhit)->hit())->cluster()->size() );
00672 //       if (dynamic_cast<const SiStripRecHit2D*>((*rhit)->hit()))      
00673 //      hChi2_vs_clsize->Fill( chi2increment, ((const SiStripRecHit2D*)(*rhit)->hit())->cluster()->amplitudes().size() );
00674       
00675 //       LogTrace("TestHits") << "rs" << std::endl;
00676 
00677 //       title.str("");
00678 //       title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
00679 //       hPullGP_X_rs[title.str()]->Fill( pullGPX_rs );
00680 //       title.str("");
00681 //       title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
00682 //       hPullGP_Y_rs[title.str()]->Fill( pullGPY_rs );
00683 //       title.str("");
00684 //       title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
00685 //       hPullGP_Z_rs[title.str()]->Fill( pullGPZ_rs );
00686 
00687 //       double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx());
00688 //       double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy());
00689 //       double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/sqrt(tsosGPEr.czz()+rhitGPEr.czz());
00690 // //       double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/*/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx())*/;
00691 // //       double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/*/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy())*/;
00692 // //       double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/*/sqrt(tsosGPEr.czz()+rhitGPEr.czz())*/;
00693 
00694 //       LogTrace("TestHits") << "tr" << std::endl;
00695 
00696 //       title.str("");
00697 //       title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
00698 //       hPullGP_X_tr[title.str()]->Fill( pullGPX_tr );
00699 //       title.str("");
00700 //       title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
00701 //       hPullGP_Y_tr[title.str()]->Fill( pullGPY_tr );
00702 //       title.str("");
00703 //       title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
00704 //       hPullGP_Z_tr[title.str()]->Fill( pullGPZ_tr );
00705 
00706 //       double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/sqrt(tsosGPEr.cxx());
00707 //       double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/sqrt(tsosGPEr.cyy());
00708 //       double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/sqrt(tsosGPEr.czz());
00709 // //       double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/*/sqrt(tsosGPEr.cxx())*/;
00710 // //       double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/*/sqrt(tsosGPEr.cyy())*/;
00711 // //       double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/*/sqrt(tsosGPEr.czz())*/;
00712 
00713 //       LogTrace("TestHits") << "ts1" << std::endl;
00714 
00715 //       title.str("");
00716 //       title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
00717 //       hPullGP_X_ts[title.str()]->Fill( pullGPX_ts );
00718 //       title.str("");
00719 //       title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
00720 //       hPullGP_Y_ts[title.str()]->Fill( pullGPY_ts );
00721 //       title.str("");
00722 //       title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
00723 //       hPullGP_Z_ts[title.str()]->Fill( pullGPZ_ts );
00724 
00725 //       double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/sqrt(tsosGMEr.cxx());
00726 //       double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/sqrt(tsosGMEr.cyy());
00727 //       double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/sqrt(tsosGMEr.czz());
00728 // //       double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/*/sqrt(tsosGMEr.cxx())*/;
00729 // //       double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/*/sqrt(tsosGMEr.cyy())*/;
00730 // //       double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/*/sqrt(tsosGMEr.czz())*/;
00731 
00732 //       LogTrace("TestHits") << "ts2" << std::endl;
00733 
00734 //       title.str("");
00735 //       title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
00736 //       hPullGM_X_ts[title.str()]->Fill( pullGMX_ts );
00737 //       title.str("");
00738 //       title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
00739 //       hPullGM_Y_ts[title.str()]->Fill( pullGMY_ts );
00740 //       title.str("");
00741 //       title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
00742 //       hPullGM_Z_ts[title.str()]->Fill( pullGMZ_ts );
00743 
00744 //       if (dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())) {
00745 //      //mono
00746 //      LogTrace("TestHits") << "MONO HIT" << std::endl;
00747 //      CTTRHp tMonoHit = 
00748 //        theBuilder->build(dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())->monoHit());
00749 //      if (tMonoHit==0) continue;
00750 //      vector<PSimHit> assMonoSimHits = hitAssociator->associateHit(*tMonoHit->hit());
00751 //      if (assMonoSimHits.size()==0) continue;
00752 //      const PSimHit sMonoHit = *(assSimHits.begin());
00753 //      const Surface * monoSurf = &( tMonoHit->det()->surface() );
00754 //      if (monoSurf==0) continue;
00755 //      TSOS monoState = thePropagator->propagate(lastState,*monoSurf);
00756 //      if (monoState.isValid()==0) continue;
00757 
00758 //      LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
00759 //      GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
00760 //      LocalPoint monoShitLPos = sMonoHit.localPosition();
00761 //      GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
00762 
00763 //      GlobalVector monoTsosGMom = monoState.globalMomentum();
00764 //      GlobalError  monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00765 //      GlobalPoint  monoTsosGPos = monoState.globalPosition();
00766 //      GlobalError  monoTsosGPEr = monoState.cartesianError().position();
00767 
00768 //      GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
00769 //      GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
00770 
00771 //      double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/sqrt(monoRhitGPEr.cxx());
00772 //      double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/sqrt(monoRhitGPEr.cyy());
00773 //      double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/sqrt(monoRhitGPEr.czz());
00774 // //   double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/*/sqrt(monoRhitGPEr.cxx())*/;
00775 // //   double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/*/sqrt(monoRhitGPEr.cyy())*/;
00776 // //   double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/*/sqrt(monoRhitGPEr.czz())*/;
00777 
00778 //      title.str("");
00779 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
00780 //      hPullGP_X_rs_mono[title.str()]->Fill( pullGPX_rs_mono );
00781 //      title.str("");
00782 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
00783 //      hPullGP_Y_rs_mono[title.str()]->Fill( pullGPY_rs_mono );
00784 //      title.str("");
00785 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
00786 //      hPullGP_Z_rs_mono[title.str()]->Fill( pullGPZ_rs_mono );
00787 
00788 //      double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx());
00789 //      double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy());
00790 //      double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz());
00791 // //   double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/*/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx())*/;
00792 // //   double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/*/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy())*/;
00793 // //   double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/*/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz())*/;
00794 
00795 //      title.str("");
00796 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
00797 //      hPullGP_X_tr_mono[title.str()]->Fill( pullGPX_tr_mono );
00798 //      title.str("");
00799 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
00800 //      hPullGP_Y_tr_mono[title.str()]->Fill( pullGPY_tr_mono );
00801 //      title.str("");
00802 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
00803 //      hPullGP_Z_tr_mono[title.str()]->Fill( pullGPZ_tr_mono );
00804 
00805 //      double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/sqrt(monoTsosGPEr.cxx());
00806 //      double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/sqrt(monoTsosGPEr.cyy());
00807 //      double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/sqrt(monoTsosGPEr.czz());
00808 // //   double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/*/sqrt(monoTsosGPEr.cxx())*/;
00809 // //   double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/*/sqrt(monoTsosGPEr.cyy())*/;
00810 // //   double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/*/sqrt(monoTsosGPEr.czz())*/;
00811 
00812 //      title.str("");
00813 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
00814 //      hPullGP_X_ts_mono[title.str()]->Fill( pullGPX_ts_mono );
00815 //      title.str("");
00816 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
00817 //      hPullGP_Y_ts_mono[title.str()]->Fill( pullGPY_ts_mono );
00818 //      title.str("");
00819 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
00820 //      hPullGP_Z_ts_mono[title.str()]->Fill( pullGPZ_ts_mono );
00821 
00822 //      double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/sqrt(monoTsosGMEr.cxx());
00823 //      double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/sqrt(monoTsosGMEr.cyy());
00824 //      double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/sqrt(monoTsosGMEr.czz());
00825 // //   double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/*/sqrt(monoTsosGMEr.cxx())*/;
00826 // //   double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/*/sqrt(monoTsosGMEr.cyy())*/;
00827 // //   double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/*/sqrt(monoTsosGMEr.czz())*/;
00828 
00829 //      title.str("");
00830 //      title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
00831 //      hPullGM_X_ts_mono[title.str()]->Fill( pullGMX_ts_mono );
00832 //      title.str("");
00833 //      title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
00834 //      hPullGM_Y_ts_mono[title.str()]->Fill( pullGMY_ts_mono );
00835 //      title.str("");
00836 //      title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
00837 //      hPullGM_Z_ts_mono[title.str()]->Fill( pullGMZ_ts_mono );
00838 
00839 //      //stereo
00840 //      LogTrace("TestHits") << "STEREO HIT" << std::endl;
00841 //      CTTRHp tStereoHit = 
00842 //        theBuilder->build(dynamic_cast<const SiStripMatchedRecHit2D*>((*rhit)->hit())->stereoHit());
00843 //      if (tStereoHit==0) continue;
00844 //      vector<PSimHit> assStereoSimHits = hitAssociator->associateHit(*tStereoHit->hit());
00845 //      if (assStereoSimHits.size()==0) continue;
00846 //      const PSimHit sStereoHit = *(assSimHits.begin());
00847 //      const Surface * stereoSurf = &( tStereoHit->det()->surface() );
00848 //      if (stereoSurf==0) continue;
00849 //      TSOS stereoState = thePropagator->propagate(lastState,*stereoSurf);
00850 //      if (stereoState.isValid()==0) continue;
00851 
00852 //      LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
00853 //      GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
00854 //      LocalPoint stereoShitLPos = sStereoHit.localPosition();
00855 //      GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
00856 
00857 //      GlobalVector stereoTsosGMom = stereoState.globalMomentum();
00858 //      GlobalError  stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00859 //      GlobalPoint  stereoTsosGPos = stereoState.globalPosition();
00860 //      GlobalError  stereoTsosGPEr = stereoState.cartesianError().position();
00861 
00862 //      GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
00863 //      GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
00864 
00865 //      double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/sqrt(stereoRhitGPEr.cxx());
00866 //      double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/sqrt(stereoRhitGPEr.cyy());
00867 //      double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/sqrt(stereoRhitGPEr.czz());
00868 // //   double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/*/sqrt(stereoRhitGPEr.cxx())*/;
00869 // //   double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/*/sqrt(stereoRhitGPEr.cyy())*/;
00870 // //   double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/*/sqrt(stereoRhitGPEr.czz())*/;
00871 
00872 //      title.str("");
00873 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
00874 //      hPullGP_X_rs_stereo[title.str()]->Fill( pullGPX_rs_stereo );
00875 //      title.str("");
00876 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
00877 //      hPullGP_Y_rs_stereo[title.str()]->Fill( pullGPY_rs_stereo );
00878 //      title.str("");
00879 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
00880 //      hPullGP_Z_rs_stereo[title.str()]->Fill( pullGPZ_rs_stereo );
00881 
00882 //      double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx());
00883 //      double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy());
00884 //      double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz());
00885 // //   double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/*/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx())*/;
00886 // //   double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/*/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy())*/;
00887 // //   double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/*/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz())*/;
00888 
00889 //      title.str("");
00890 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
00891 //      hPullGP_X_tr_stereo[title.str()]->Fill( pullGPX_tr_stereo );
00892 //      title.str("");
00893 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
00894 //      hPullGP_Y_tr_stereo[title.str()]->Fill( pullGPY_tr_stereo );
00895 //      title.str("");
00896 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
00897 //      hPullGP_Z_tr_stereo[title.str()]->Fill( pullGPZ_tr_stereo );
00898 
00899 //      double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/sqrt(stereoTsosGPEr.cxx());
00900 //      double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/sqrt(stereoTsosGPEr.cyy());
00901 //      double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/sqrt(stereoTsosGPEr.czz());
00902 // //   double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/*/sqrt(stereoTsosGPEr.cxx())*/;
00903 // //   double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/*/sqrt(stereoTsosGPEr.cyy())*/;
00904 // //   double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/*/sqrt(stereoTsosGPEr.czz())*/;
00905 
00906 //      title.str("");
00907 //      title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
00908 //      hPullGP_X_ts_stereo[title.str()]->Fill( pullGPX_ts_stereo );
00909 //      title.str("");
00910 //      title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00911 //      hPullGP_Y_ts_stereo[title.str()]->Fill( pullGPY_ts_stereo );
00912 //      title.str("");
00913 //      title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00914 //      hPullGP_Z_ts_stereo[title.str()]->Fill( pullGPZ_ts_stereo );
00915 
00916 //      double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/sqrt(stereoTsosGMEr.cxx());
00917 //      double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/sqrt(stereoTsosGMEr.cyy());
00918 //      double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/sqrt(stereoTsosGMEr.czz());
00919 // //   double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/*/sqrt(stereoTsosGMEr.cxx())*/;
00920 // //   double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/*/sqrt(stereoTsosGMEr.cyy())*/;
00921 // //   double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/*/sqrt(stereoTsosGMEr.czz())*/;
00922 
00923 //      title.str("");
00924 //      title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
00925 //      hPullGM_X_ts_stereo[title.str()]->Fill( pullGMX_ts_stereo );
00926 //      title.str("");
00927 //      title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00928 //      hPullGM_Y_ts_stereo[title.str()]->Fill( pullGMY_ts_stereo );
00929 //      title.str("");
00930 //      title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00931 //      hPullGM_Z_ts_stereo[title.str()]->Fill( pullGMZ_ts_stereo );
00932 //       }
00933 //     }
00934 //   }
00935 //   delete hitAssociator;
00936 //   LogTrace("TestHits") << "end of event" << std::endl;
00937 // }
00938 
00939 void TestHits::endJob() {
00940   //file->Write();
00941   TDirectory * chi2i = file->mkdir("Chi2_Increment");
00942 
00943   TDirectory * gp_ts = file->mkdir("GP_TSOS-SimHit");
00944   TDirectory * gm_ts = file->mkdir("GM_TSOS-SimHit");
00945   TDirectory * gp_tr = file->mkdir("GP_TSOS-RecHit");
00946   TDirectory * gp_rs = file->mkdir("GP_RecHit-SimHit");
00947 
00948   TDirectory * gp_tsx = gp_ts->mkdir("X");
00949   TDirectory * gp_tsy = gp_ts->mkdir("Y");
00950   TDirectory * gp_tsz = gp_ts->mkdir("Z");
00951   TDirectory * gm_tsx = gm_ts->mkdir("X");
00952   TDirectory * gm_tsy = gm_ts->mkdir("Y");
00953   TDirectory * gm_tsz = gm_ts->mkdir("Z");
00954   TDirectory * gp_trx = gp_tr->mkdir("X");
00955   TDirectory * gp_try = gp_tr->mkdir("Y");
00956   TDirectory * gp_trz = gp_tr->mkdir("Z");
00957   TDirectory * gp_rsx = gp_rs->mkdir("X");
00958   TDirectory * gp_rsy = gp_rs->mkdir("Y");
00959   TDirectory * gp_rsz = gp_rs->mkdir("Z");
00960 
00961   TDirectory * gp_tsx_mono = gp_ts->mkdir("MONOX");
00962   TDirectory * gp_tsy_mono = gp_ts->mkdir("MONOY");
00963   TDirectory * gp_tsz_mono = gp_ts->mkdir("MONOZ");
00964   TDirectory * gm_tsx_mono = gm_ts->mkdir("MONOX");
00965   TDirectory * gm_tsy_mono = gm_ts->mkdir("MONOY");
00966   TDirectory * gm_tsz_mono = gm_ts->mkdir("MONOZ");
00967   TDirectory * gp_trx_mono = gp_tr->mkdir("MONOX");
00968   TDirectory * gp_try_mono = gp_tr->mkdir("MONOY");
00969   TDirectory * gp_trz_mono = gp_tr->mkdir("MONOZ");
00970   TDirectory * gp_rsx_mono = gp_rs->mkdir("MONOX");
00971   TDirectory * gp_rsy_mono = gp_rs->mkdir("MONOY");
00972   TDirectory * gp_rsz_mono = gp_rs->mkdir("MONOZ");
00973 
00974   TDirectory * gp_tsx_stereo = gp_ts->mkdir("STEREOX");
00975   TDirectory * gp_tsy_stereo = gp_ts->mkdir("STEREOY");
00976   TDirectory * gp_tsz_stereo = gp_ts->mkdir("STEREOZ");
00977   TDirectory * gm_tsx_stereo = gm_ts->mkdir("STEREOX");
00978   TDirectory * gm_tsy_stereo = gm_ts->mkdir("STEREOY");
00979   TDirectory * gm_tsz_stereo = gm_ts->mkdir("STEREOZ");
00980   TDirectory * gp_trx_stereo = gp_tr->mkdir("STEREOX");
00981   TDirectory * gp_try_stereo = gp_tr->mkdir("STEREOY");
00982   TDirectory * gp_trz_stereo = gp_tr->mkdir("STEREOZ");
00983   TDirectory * gp_rsx_stereo = gp_rs->mkdir("STEREOX");
00984   TDirectory * gp_rsy_stereo = gp_rs->mkdir("STEREOY");
00985   TDirectory * gp_rsz_stereo = gp_rs->mkdir("STEREOZ");
00986 
00987   chi2i->cd();
00988   hTotChi2Increment->Write();
00989   hProcess_vs_Chi2->Write();
00990   hClsize_vs_Chi2->Write();
00991   for (int i=0; i!=6; i++)
00992     for (int j=0; j!=9; j++){
00993       if (i==0 && j>2) break;
00994       if (i==1 && j>1) break;
00995       if (i==2 && j>3) break;
00996       if (i==3 && j>2) break;
00997       if (i==4 && j>5) break;
00998       if (i==5 && j>8) break;
00999       chi2i->cd();
01000       title.str("");
01001       title << "Chi2Increment_" << i+1 << "-" << j+1;
01002       hChi2Increment[title.str()]->Write();
01003 
01004       gp_ts->cd();
01005       gp_tsx->cd();
01006       title.str("");
01007       title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
01008       hPullGP_X_ts[title.str()]->Write();
01009       gp_tsy->cd();
01010       title.str("");
01011       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts";
01012       hPullGP_Y_ts[title.str()]->Write();
01013       gp_tsz->cd();
01014       title.str("");
01015       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts";
01016       hPullGP_Z_ts[title.str()]->Write();
01017 
01018       gm_ts->cd();
01019       gm_tsx->cd();
01020       title.str("");
01021       title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
01022       hPullGM_X_ts[title.str()]->Write();
01023       gm_tsy->cd();
01024       title.str("");
01025       title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
01026       hPullGM_Y_ts[title.str()]->Write();
01027       gm_tsz->cd();
01028       title.str("");
01029       title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
01030       hPullGM_Z_ts[title.str()]->Write();
01031 
01032       gp_tr->cd();
01033       gp_trx->cd();
01034       title.str("");
01035       title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
01036       hPullGP_X_tr[title.str()]->Write();
01037       gp_try->cd();
01038       title.str("");
01039       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
01040       hPullGP_Y_tr[title.str()]->Write();
01041       gp_trz->cd();
01042       title.str("");
01043       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
01044       hPullGP_Z_tr[title.str()]->Write();
01045 
01046       gp_rs->cd();
01047       gp_rsx->cd();
01048       title.str("");
01049       title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
01050       hPullGP_X_rs[title.str()]->Write();
01051       gp_rsy->cd();
01052       title.str("");
01053       title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
01054       hPullGP_Y_rs[title.str()]->Write();
01055       gp_rsz->cd();
01056       title.str("");
01057       title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
01058       hPullGP_Z_rs[title.str()]->Write();
01059 
01060       if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
01061         //mono
01062         gp_ts->cd();
01063         gp_tsx_mono->cd();
01064         title.str("");
01065         title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
01066         hPullGP_X_ts_mono[title.str()]->Write();
01067         gp_tsy_mono->cd();
01068         title.str("");
01069         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
01070         hPullGP_Y_ts_mono[title.str()]->Write();
01071         gp_tsz_mono->cd();
01072         title.str("");
01073         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
01074         hPullGP_Z_ts_mono[title.str()]->Write();
01075 
01076         gm_ts->cd();
01077         gm_tsx_mono->cd();
01078         title.str("");
01079         title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
01080         hPullGM_X_ts_mono[title.str()]->Write();
01081         gm_tsy_mono->cd();
01082         title.str("");
01083         title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
01084         hPullGM_Y_ts_mono[title.str()]->Write();
01085         gm_tsz_mono->cd();
01086         title.str("");
01087         title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
01088         hPullGM_Z_ts_mono[title.str()]->Write();
01089 
01090         gp_tr->cd();
01091         gp_trx_mono->cd();
01092         title.str("");
01093         title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
01094         hPullGP_X_tr_mono[title.str()]->Write();
01095         gp_try_mono->cd();
01096         title.str("");
01097         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
01098         hPullGP_Y_tr_mono[title.str()]->Write();
01099         gp_trz_mono->cd();
01100         title.str("");
01101         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
01102         hPullGP_Z_tr_mono[title.str()]->Write();
01103 
01104         gp_rs->cd();
01105         gp_rsx_mono->cd();
01106         title.str("");
01107         title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
01108         hPullGP_X_rs_mono[title.str()]->Write();
01109         gp_rsy_mono->cd();
01110         title.str("");
01111         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
01112         hPullGP_Y_rs_mono[title.str()]->Write();
01113         gp_rsz_mono->cd();
01114         title.str("");
01115         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
01116         hPullGP_Z_rs_mono[title.str()]->Write();
01117 
01118         //stereo
01119         gp_ts->cd();
01120         gp_tsx_stereo->cd();
01121         title.str("");
01122         title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
01123         hPullGP_X_ts_stereo[title.str()]->Write();
01124         gp_tsy_stereo->cd();
01125         title.str("");
01126         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
01127         hPullGP_Y_ts_stereo[title.str()]->Write();
01128         gp_tsz_stereo->cd();
01129         title.str("");
01130         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
01131         hPullGP_Z_ts_stereo[title.str()]->Write();
01132 
01133         gm_ts->cd();
01134         gm_tsx_stereo->cd();
01135         title.str("");
01136         title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
01137         hPullGM_X_ts_stereo[title.str()]->Write();
01138         gm_tsy_stereo->cd();
01139         title.str("");
01140         title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
01141         hPullGM_Y_ts_stereo[title.str()]->Write();
01142         gm_tsz_stereo->cd();
01143         title.str("");
01144         title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
01145         hPullGM_Z_ts_stereo[title.str()]->Write();
01146 
01147         gp_tr->cd();
01148         gp_trx_stereo->cd();
01149         title.str("");
01150         title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
01151         hPullGP_X_tr_stereo[title.str()]->Write();
01152         gp_try_stereo->cd();
01153         title.str("");
01154         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
01155         hPullGP_Y_tr_stereo[title.str()]->Write();
01156         gp_trz_stereo->cd();
01157         title.str("");
01158         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
01159         hPullGP_Z_tr_stereo[title.str()]->Write();
01160 
01161         gp_rs->cd();
01162         gp_rsx_stereo->cd();
01163         title.str("");
01164         title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
01165         hPullGP_X_rs_stereo[title.str()]->Write();
01166         gp_rsy_stereo->cd();
01167         title.str("");
01168         title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
01169         hPullGP_Y_rs_stereo[title.str()]->Write();
01170         gp_rsz_stereo->cd();
01171         title.str("");
01172         title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
01173         hPullGP_Z_rs_stereo[title.str()]->Write();
01174       }
01175     }
01176 
01177   file->Close();
01178 }
01179 
01180 //needed by to do the residual for matched hits
01181 //taken from SiStripTrackingRecHitsValid.cc
01182 std::pair<LocalPoint,LocalVector> 
01183 TestHits::projectHit( const PSimHit& hit, const StripGeomDetUnit* stripDet, const BoundPlane& plane) 
01184 { 
01185    const StripTopology& topol = stripDet->specificTopology();
01186    GlobalPoint globalpos= stripDet->surface().toGlobal(hit.localPosition());
01187    LocalPoint localHit = plane.toLocal(globalpos);
01188    //track direction
01189    LocalVector locdir=hit.localDirection();
01190    //rotate track in new frame
01191    
01192    GlobalVector globaldir= stripDet->surface().toGlobal(locdir);
01193    LocalVector dir=plane.toLocal(globaldir);
01194    float scale = -localHit.z() / dir.z();
01195    
01196    LocalPoint projectedPos = localHit + scale*dir;
01197    
01198    float selfAngle = topol.stripAngle( topol.strip( hit.localPosition()));
01199    
01200    LocalVector stripDir( sin(selfAngle), cos(selfAngle), 0); // vector along strip in hit frame
01201    
01202    LocalVector localStripDir( plane.toLocal(stripDet->surface().toGlobal( stripDir)));
01203    
01204    return std::pair<LocalPoint,LocalVector>( projectedPos, localStripDir);
01205 }
01206 
01207 #include "FWCore/Framework/interface/ModuleFactory.h"
01208 #include "FWCore/Framework/interface/MakerMacros.h"
01209 
01210 DEFINE_FWK_MODULE(TestHits);
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_7/src/RecoTracker/DebugTools/plugins/TestHits.cc
