CMS 3D CMS Logo
RecoTracker » DebugTools » plugins

TestHits.cc

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