CMS 3D CMS Logo

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