00001 #include "RecoTracker/DebugTools/interface/TestTrackHits.h"
00002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00003 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00004 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
00005 #include "TrackingTools/MaterialEffects/interface/PropagatorWithMaterial.h"
00006 #include "TrackingTools/TrackFitters/interface/TrajectoryStateCombiner.h"
00007 #include "SimTracker/Records/interface/TrackAssociatorRecord.h"
00008 #include <TDirectory.h>
00009
00010 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
00011 #include "Geometry/Records/interface/IdealGeometryRecord.h"
00012 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00013
00014 typedef TrajectoryStateOnSurface TSOS;
00015 typedef TransientTrackingRecHit::ConstRecHitPointer CTTRHp;
00016 using namespace std;
00017 using namespace edm;
00018
00019 TestTrackHits::TestTrackHits(const edm::ParameterSet& iConfig):
00020 conf_(iConfig) {
00021 LogTrace("TestTrackHits") << conf_;
00022 propagatorName = conf_.getParameter<std::string>("Propagator");
00023 builderName = conf_.getParameter<std::string>("TTRHBuilder");
00024 srcName = conf_.getParameter<std::string>("src");
00025 tpName = conf_.getParameter<std::string>("tpname");
00026 updatorName = conf_.getParameter<std::string>("updator");
00027 out = conf_.getParameter<std::string>("out");
00028
00029
00030
00031
00032
00033
00034
00035
00036 }
00037
00038 TestTrackHits::~TestTrackHits(){}
00039
00040 void TestTrackHits::beginRun(edm::Run & run, const edm::EventSetup& iSetup)
00041 {
00042 iSetup.get<TrackerDigiGeometryRecord>().get(theG);
00043 iSetup.get<IdealMagneticFieldRecord>().get(theMF);
00044 iSetup.get<TrackingComponentsRecord>().get(propagatorName,thePropagator);
00045 iSetup.get<TransientRecHitRecord>().get(builderName,theBuilder);
00046 iSetup.get<TrackingComponentsRecord>().get(updatorName,theUpdator);
00047 iSetup.get<TrackAssociatorRecord>().get("TrackAssociatorByHits",trackAssociator);
00048
00049 file = new TFile(out.c_str(),"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
00059 title.str("");
00060 title << "Chi2Increment_" << i+1 << "-" << j+1 ;
00061 hChi2Increment[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00062 title.str("");
00063 title << "Chi2IncrementVsEta_" << i+1 << "-" << j+1 ;
00064 hChi2IncrementVsEta[title.str()] = new TH2F(title.str().c_str(),title.str().c_str(),50,-2.5,2.5,1000,0,100);
00065 title.str("");
00066 title << "Chi2GoodHit_" << i+1 << "-" << j+1 ;
00067 hChi2GoodHit[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00068 title.str("");
00069 title << "Chi2BadHit_" << i+1 << "-" << j+1 ;
00070 hChi2BadHit[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00071 title.str("");
00072 title << "Chi2DeltaHit_" << i+1 << "-" << j+1 ;
00073 hChi2DeltaHit[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00074 title.str("");
00075 title << "Chi2NSharedHit_" << i+1 << "-" << j+1 ;
00076 hChi2NSharedHit[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00077 title.str("");
00078 title << "Chi2SharedHit_" << i+1 << "-" << j+1 ;
00079 hChi2SharedHit[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00080
00081 title.str("");
00082 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
00083 hPullGP_X_ts[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 << "_ts";
00086 hPullGP_Y_ts[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 << "_ts";
00089 hPullGP_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00090
00091 title.str("");
00092 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
00093 hPullGM_X_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00094 title.str("");
00095 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
00096 hPullGM_Y_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00097 title.str("");
00098 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
00099 hPullGM_Z_ts[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00100
00101 title.str("");
00102 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
00103 hPullGP_X_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00104 title.str("");
00105 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
00106 hPullGP_Y_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00107 title.str("");
00108 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
00109 hPullGP_Z_tr[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00110
00111 title.str("");
00112 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
00113 hPullGP_X_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00114 title.str("");
00115 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
00116 hPullGP_Y_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00117 title.str("");
00118 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
00119 hPullGP_Z_rs[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00120
00121 if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
00122
00123 title.str("");
00124 title << "Chi2Increment_mono_" << i+1 << "-" << j+1 ;
00125 hChi2Increment_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00126
00127 title.str("");
00128 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
00129 hPullGP_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00130 title.str("");
00131 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
00132 hPullGP_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00133 title.str("");
00134 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
00135 hPullGP_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00136
00137 title.str("");
00138 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
00139 hPullGM_X_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00140 title.str("");
00141 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
00142 hPullGM_Y_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00143 title.str("");
00144 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
00145 hPullGM_Z_ts_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00146
00147 title.str("");
00148 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
00149 hPullGP_X_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00150 title.str("");
00151 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
00152 hPullGP_Y_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00153 title.str("");
00154 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
00155 hPullGP_Z_tr_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00156
00157 title.str("");
00158 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
00159 hPullGP_X_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00160 title.str("");
00161 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
00162 hPullGP_Y_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00163 title.str("");
00164 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
00165 hPullGP_Z_rs_mono[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00166
00167
00168 title.str("");
00169 title << "Chi2Increment_stereo_" << i+1 << "-" << j+1 ;
00170 hChi2Increment_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,0,100);
00171
00172 title.str("");
00173 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
00174 hPullGP_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00175 title.str("");
00176 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
00177 hPullGP_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00178 title.str("");
00179 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
00180 hPullGP_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00181
00182 title.str("");
00183 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
00184 hPullGM_X_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00185 title.str("");
00186 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
00187 hPullGM_Y_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00188 title.str("");
00189 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
00190 hPullGM_Z_ts_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00191
00192 title.str("");
00193 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
00194 hPullGP_X_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00195 title.str("");
00196 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
00197 hPullGP_Y_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00198 title.str("");
00199 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
00200 hPullGP_Z_tr_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00201
00202 title.str("");
00203 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
00204 hPullGP_X_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00205 title.str("");
00206 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
00207 hPullGP_Y_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00208 title.str("");
00209 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
00210 hPullGP_Z_rs_stereo[title.str()] = new TH1F(title.str().c_str(),title.str().c_str(),1000,-50,50);
00211 }
00212 }
00213 hTotChi2Increment = new TH1F("TotChi2Increment","TotChi2Increment",1000,0,100);
00214 hTotChi2GoodHit = new TH1F("TotChi2GoodHit","TotChi2GoodHit",1000,0,100);
00215 hTotChi2BadHit = new TH1F("TotChi2BadHit","TotChi2BadHit",1000,0,100);
00216 hTotChi2DeltaHit = new TH1F("TotChi2DeltaHit","TotChi2DeltaHit",1000,0,100);
00217 hTotChi2NSharedHit = new TH1F("TotChi2NSharedHit","TotChi2NSharedHit",1000,0,100);
00218 hTotChi2SharedHit = new TH1F("TotChi2SharedHit","TotChi2SharedHit",1000,0,100);
00219 hProcess_vs_Chi2 = new TH2F("Process_vs_Chi2","Process_vs_Chi2",1000,0,100,17,-0.5,16.5);
00220 hClsize_vs_Chi2 = new TH2F("Clsize_vs_Chi2","Clsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00221 hPixClsize_vs_Chi2= new TH2F("PixClsize_vs_Chi2","PixClsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00222 hPrjClsize_vs_Chi2= new TH2F("PrjClsize_vs_Chi2","PrjClsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00223 hSt1Clsize_vs_Chi2= new TH2F("St1Clsize_vs_Chi2","St1Clsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00224 hSt2Clsize_vs_Chi2= new TH2F("St2Clsize_vs_Chi2","St2Clsize_vs_Chi2",1000,0,100,17,-0.5,16.5);
00225 hGoodHit_vs_Chi2 = new TH2F("GoodHit_vs_Chi2","GoodHit_vs_Chi2",10000,0,1000,2,-0.5,1.5);
00226 hClusterSize = new TH1F("ClusterSize","ClusterSize",40,-0.5,39.5);
00227 hPixClusterSize = new TH1F("PixClusterSize","PixClusterSize",40,-0.5,39.5);
00228 hPrjClusterSize = new TH1F("PrjClusterSize","PrjClusterSize",40,-0.5,39.5);
00229 hSt1ClusterSize = new TH1F("St1ClusterSize","St1ClusterSize",40,-0.5,39.5);
00230 hSt2ClusterSize = new TH1F("St2ClusterSize","St2ClusterSize",40,-0.5,39.5);
00231 hSimHitVecSize = new TH1F("hSimHitVecSize","hSimHitVecSize",40,-0.5,39.5);
00232 hPixSimHitVecSize = new TH1F("PixSimHitVecSize","PixSimHitVecSize",40,-0.5,39.5);
00233 hPrjSimHitVecSize = new TH1F("PrjSimHitVecSize","PrjSimHitVecSize",40,-0.5,39.5);
00234 hSt1SimHitVecSize = new TH1F("St1SimHitVecSize","St1SimHitVecSize",40,-0.5,39.5);
00235 hSt2SimHitVecSize = new TH1F("St2SimHitVecSize","St2SimHitVecSize",40,-0.5,39.5);
00236 goodbadmerged = new TH1F("goodbadmerged","goodbadmerged",5,0.5,5.5);
00237 energyLossRatio = new TH1F("energyLossRatio","energyLossRatio",100,0,1);
00238 mergedPull = new TH1F("mergedPull","mergedPull",200,0,20);
00239 probXgood = new TH1F("probXgood","probXgood",110,0,1.1);
00240 probXbad = new TH1F("probXbad","probXbad",110,0,1.1);
00241 probXdelta = new TH1F("probXdelta","probXdelta",110,0,1.1);
00242 probXshared = new TH1F("probXshared","probXshared",110,0,1.1);
00243 probXnoshare = new TH1F("probXnoshare","probXnoshare",110,0,1.1);
00244 probYgood = new TH1F("probYgood","probYgood",110,0,1.1);
00245 probYbad = new TH1F("probYbad","probYbad",110,0,1.1);
00246 probYdelta = new TH1F("probYdelta","probYdelta",110,0,1.1);
00247 probYshared = new TH1F("probYshared","probYshared",110,0,1.1);
00248 probYnoshare = new TH1F("probYnoshare","probYnoshare",110,0,1.1);
00249 }
00250
00251 void TestTrackHits::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00252 {
00253
00254 edm::ESHandle<TrackerTopology> tTopo;
00255 iSetup.get<IdealGeometryRecord>().get(tTopo);
00256
00257
00258 LogDebug("TestTrackHits") << "new event" ;
00259 iEvent.getByLabel(srcName,trajCollectionHandle);
00260 iEvent.getByLabel(srcName,trackCollectionHandle);
00261 iEvent.getByLabel(srcName,trajTrackAssociationCollectionHandle);
00262 iEvent.getByLabel(tpName,trackingParticleCollectionHandle);
00263 edm::Handle<reco::BeamSpot> beamSpot;
00264 iEvent.getByLabel("offlineBeamSpot",beamSpot);
00265
00266 LogTrace("TestTrackHits") << "Tr collection size=" << trackCollectionHandle->size();
00267 LogTrace("TestTrackHits") << "TP collection size=" << trackingParticleCollectionHandle->size();
00268
00269 hitAssociator = new TrackerHitAssociator(iEvent);
00270
00271 reco::RecoToSimCollection recSimColl=trackAssociator->associateRecoToSim(trackCollectionHandle,
00272 trackingParticleCollectionHandle,
00273 &iEvent,&iSetup);
00274
00275 TrajectoryStateCombiner combiner;
00276
00277 int evtHits = 0;
00278 int i=0;
00279 int yy=0;
00280 int yyy=0;
00281 for(std::vector<Trajectory>::const_iterator it = trajCollectionHandle->begin(); it!=trajCollectionHandle->end();it++){
00282
00283 LogTrace("TestTrackHits") << "\n*****************new trajectory********************" ;
00284 double tchi2 = 0;
00285
00286 std::vector<TrajectoryMeasurement> tmColl = it->measurements();
00287
00288 edm::Ref<std::vector<Trajectory> > traj(trajCollectionHandle, i);
00289 reco::TrackRef tmptrack = (*trajTrackAssociationCollectionHandle.product())[traj];
00290 edm::RefToBase<reco::Track> track(tmptrack);
00291
00292
00293
00294
00295
00296
00297 std::vector<std::pair<TrackingParticleRef, double> > tP;
00298 if(recSimColl.find(track) != recSimColl.end()){
00299 tP = recSimColl[track];
00300 if (tP.size()!=0) {
00301 edm::LogVerbatim("TestTrackHits") << "reco::Track #" << ++yyy << " with pt=" << track->pt()
00302 << " associated with quality:" << tP.begin()->second <<" good track #" << ++yy << " has hits:" << track->numberOfValidHits() << "\n";
00303 }
00304 }else{
00305 edm::LogVerbatim("TestTrackHits") << "reco::Track #" << ++yyy << " with pt=" << track->pt()
00306 << " NOT associated to any TrackingParticle" << "\n";
00307 i++;
00308 continue;
00309 }
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320 TrackingParticleRef tp = tP.begin()->first;
00321 LogTrace("TestTrackHits") << "a tp is associated with fraction=" << tP.begin()->second;
00322
00323 std::vector<unsigned int> tpids;
00324 for (TrackingParticle::g4t_iterator g4T=tp->g4Track_begin(); g4T!=tp->g4Track_end(); ++g4T) {
00325 LogTrace("TestTrackHits") << "tp id=" << g4T->trackId();
00326 tpids.push_back(g4T->trackId());
00327 }
00328
00329
00330 int pp = 0;
00331 for (std::vector<TrajectoryMeasurement>::iterator tm=tmColl.begin();tm!=tmColl.end();++tm){
00332
00333 tchi2+=tm->estimate();
00334
00335 LogTrace("TestTrackHits") << "+++++++++++++++++new hit+++++++++++++++++" ;
00336 CTTRHp rhit = tm->recHit();
00337
00338
00339 TSOS state = combiner(tm->backwardPredictedState(), tm->forwardPredictedState());
00340
00341 if (rhit->isValid()==0 && rhit->det()!=0) continue;
00342 evtHits++;
00343 LogTrace("TestTrackHits") << "valid hit #" << ++pp << "of hits=" << track->numberOfValidHits();
00344
00345 int subdetId = rhit->det()->geographicalId().subdetId();
00346 DetId id = rhit->det()->geographicalId();
00347 int layerId = tTopo->layer(id);
00348 LogTrace("TestTrackHits") << "subdetId=" << subdetId << " layerId=" << layerId ;
00349
00350 const Surface * surf = rhit->surface();
00351 if (surf==0) continue;
00352
00353 double energyLoss_ = 0.;
00354 unsigned int monoId = 0;
00355 std::vector<double> energyLossM;
00356 std::vector<double> energyLossS;
00357 std::vector<PSimHit> assSimHits = hitAssociator->associateHit(*(rhit)->hit());
00358 unsigned int simhitvecsize = assSimHits.size();
00359 if (simhitvecsize==0) continue;
00360 PSimHit shit;
00361 std::vector<unsigned int> trackIds;
00362 energyLossS.clear();
00363 energyLossM.clear();
00364 for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
00365 unsigned int tId = m->trackId();
00366 if (find(trackIds.begin(),trackIds.end(),tId)==trackIds.end()) trackIds.push_back(tId);
00367 if (m->energyLoss()>energyLoss_) {
00368 shit=*m;
00369 energyLoss_ = m->energyLoss();
00370 }
00371 if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())) {
00372 if (monoId==0) monoId = m->detUnitId();
00373 if (monoId == m->detUnitId()){
00374 energyLossM.push_back(m->energyLoss());
00375 }
00376 else {
00377 energyLossS.push_back(m->energyLoss());
00378 }
00379
00380 } else {
00381 energyLossM.push_back(m->energyLoss());
00382 }
00383 }
00384
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398 double chi2increment = tm->estimate();
00399 LogTrace("TestTrackHits") << "tm->estimate()=" << tm->estimate();
00400 title.str("");
00401 title << "Chi2Increment_" << subdetId << "-" << layerId;
00402 hChi2Increment[title.str()]->Fill( chi2increment );
00403 title.str("");
00404 title << "Chi2IncrementVsEta_" << subdetId << "-" << layerId;
00405 hChi2IncrementVsEta[title.str()]->Fill( track->eta(), chi2increment );
00406 hTotChi2Increment->Fill( chi2increment );
00407 hProcess_vs_Chi2->Fill( chi2increment, shit.processType() );
00408
00409 int clustersize = 0;
00410 bool mergedhit = false;
00411 if (dynamic_cast<const SiPixelRecHit*>(rhit->hit())){
00412 clustersize = ((const SiPixelRecHit*)(rhit->hit()))->cluster()->size() ;
00413 hPixClsize_vs_Chi2->Fill(chi2increment, clustersize);
00414 hPixClusterSize->Fill(clustersize);
00415 hPixSimHitVecSize->Fill(simhitvecsize);
00416 if (simhitvecsize>1) mergedhit = true;
00417 }
00418 else if (dynamic_cast<const SiStripRecHit2D*>(rhit->hit())){
00419 clustersize = ((const SiStripRecHit2D*)(rhit->hit()))->cluster()->amplitudes().size() ;
00420 hSt1Clsize_vs_Chi2->Fill(chi2increment, clustersize);
00421 hSt1ClusterSize->Fill(clustersize);
00422 hSt1SimHitVecSize->Fill(simhitvecsize);
00423 if (simhitvecsize>1) mergedhit = true;
00424 }
00425 else if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())){
00426 int clsize1 = ((const SiStripMatchedRecHit2D*)(rhit->hit()))->monoCluster().amplitudes().size();
00427 int clsize2 = ((const SiStripMatchedRecHit2D*)(rhit->hit()))->stereoCluster().amplitudes().size();
00428 if (clsize1>clsize2) clustersize = clsize1;
00429 else clustersize = clsize2;
00430 hSt2Clsize_vs_Chi2->Fill(chi2increment, clustersize);
00431 hSt2ClusterSize->Fill(clustersize);
00432 hSt2SimHitVecSize->Fill(simhitvecsize);
00433 if (simhitvecsize>2) mergedhit = true;
00434 }
00435 else if (dynamic_cast<const ProjectedSiStripRecHit2D*>(rhit->hit())){
00436 clustersize = ((const ProjectedSiStripRecHit2D*)(rhit->hit()))->originalHit().cluster()->amplitudes().size();
00437 hPrjClsize_vs_Chi2->Fill(chi2increment, clustersize);
00438 hPrjClusterSize->Fill(clustersize);
00439 hPrjSimHitVecSize->Fill(simhitvecsize);
00440 if (simhitvecsize>1) mergedhit = true;
00441 }
00442 hClsize_vs_Chi2->Fill( chi2increment, clustersize);
00443 hClusterSize->Fill(clustersize);
00444 hSimHitVecSize->Fill(simhitvecsize);
00445
00446
00447
00448
00449
00450
00451 std::vector<SimHitIdpr> simTrackIds = hitAssociator->associateHitId(*(rhit)->hit());
00452 bool goodhit = false;
00453 for(size_t j=0; j<simTrackIds.size(); j++){
00454 LogTrace("TestTrackHits") << "hit id=" << simTrackIds[j].first;
00455 for (size_t jj=0; jj<tpids.size(); jj++){
00456 if (simTrackIds[j].first == tpids[jj]) goodhit = true;
00457 break;
00458 }
00459 }
00460 bool shared = true;
00461 bool ioniOnly = true;
00462 const SiPixelRecHit* pix = dynamic_cast<const SiPixelRecHit*>(rhit->hit());
00463 if (goodhit) {
00464 if (energyLossM.size()>1&&energyLossS.size()<=1) {
00465 sort(energyLossM.begin(),energyLossM.end(),greater<double>());
00466 energyLossRatio->Fill(energyLossM[1]/energyLossM[0]);
00467 }
00468 else if (energyLossS.size()>1&&energyLossM.size()<=1) {
00469 sort(energyLossS.begin(),energyLossS.end(),greater<double>());
00470 energyLossRatio->Fill(energyLossS[1]/energyLossS[0]);
00471 }
00472 else if (energyLossS.size()>1&&energyLossM.size()>1) {
00473 sort(energyLossM.begin(),energyLossM.end(),greater<double>());
00474 sort(energyLossS.begin(),energyLossS.end(),greater<double>());
00475 energyLossM[1]/energyLossM[0] > energyLossS[1]/energyLossS[0]
00476 ? energyLossRatio->Fill(energyLossM[1]/energyLossM[0])
00477 : energyLossRatio->Fill(energyLossS[1]/energyLossS[0]);
00478 }
00479 if ( mergedhit ) {
00480
00481 LogVerbatim("TestTrackHits") << "MERGED HIT" << std::endl;
00482 unsigned int idc = 0;
00483 for (size_t jj=0; jj<tpids.size(); jj++){
00484 idc += std::count(trackIds.begin(),trackIds.end(),tpids[jj]);
00485 }
00486 if (idc==trackIds.size()) {
00487 shared = false;
00488 }
00489 for(std::vector<PSimHit>::const_iterator m=assSimHits.begin()+1; m<assSimHits.end(); m++){
00490 if ((m->processType()!=7&&m->processType()!=8&&m->processType()!=9)&&abs(m->particleType())!=11){
00491 ioniOnly = false;
00492 break;
00493 }
00494 }
00495 if (ioniOnly&&!shared) {
00496 title.str("");
00497 title << "Chi2DeltaHit_" << subdetId << "-" << layerId;
00498 hChi2DeltaHit[title.str()]->Fill( chi2increment );
00499 hTotChi2DeltaHit->Fill( chi2increment );
00500 if (pix) {
00501 probXdelta->Fill(pix->probabilityX());
00502 probYdelta->Fill(pix->probabilityY());
00503 }
00504 } else if(!ioniOnly&&!shared) {
00505 title.str("");
00506 title << "Chi2NSharedHit_" << subdetId << "-" << layerId;
00507 hChi2NSharedHit[title.str()]->Fill( chi2increment );
00508 hTotChi2NSharedHit->Fill( chi2increment );
00509 if (pix) {
00510 probXnoshare->Fill(pix->probabilityX());
00511 probYnoshare->Fill(pix->probabilityY());
00512 }
00513 } else {
00514 title.str("");
00515 title << "Chi2SharedHit_" << subdetId << "-" << layerId;
00516 hChi2SharedHit[title.str()]->Fill( chi2increment );
00517 hTotChi2SharedHit->Fill( chi2increment );
00518 if (pix) {
00519 probXshared->Fill(pix->probabilityX());
00520 probYshared->Fill(pix->probabilityY());
00521 }
00522 }
00523
00524 for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++) {
00525 unsigned int tId = m->trackId();
00526 LogVerbatim("TestTrackHits") << "component with id=" << tId << " eLoss=" << m->energyLoss() << " pType=" << m->processType();
00527 if (find(tpids.begin(),tpids.end(),tId)==tpids.end()) continue;
00528 if (m->processType()==2) {
00529 GlobalPoint gpr = rhit->globalPosition();
00530 AlgebraicSymMatrix33 ger = rhit->globalPositionError().matrix();
00531 GlobalPoint gps = surf->toGlobal(m->localPosition());
00532 LogVerbatim("TestTrackHits") << gpr << " " << gps << " " << ger;
00533 ROOT::Math::SVector<double,3> delta;
00534 delta[0]=gpr.x()-gps.x();
00535 delta[1]=gpr.y()-gps.y();
00536 delta[2]=gpr.z()-gps.z();
00537 LogVerbatim("TestTrackHits") << delta << " " << ger ;
00538 double mpull = sqrt(delta[0]*delta[0]/ger[0][0]+delta[1]*delta[1]/ger[1][1]+delta[2]*delta[2]/ger[2][2]);
00539 LogVerbatim("TestTrackHits") << "hit pull=" << mpull;
00540 mergedPull->Fill(mpull);
00541 break;
00542 }
00543 }
00544 } else {
00545 LogVerbatim("TestTrackHits") << "good hit" ;
00546 title.str("");
00547 title << "Chi2GoodHit_" << subdetId << "-" << layerId;
00548 hChi2GoodHit[title.str()]->Fill( chi2increment );
00549 hTotChi2GoodHit->Fill( chi2increment );
00550 if (pix) {
00551 probXgood->Fill(pix->probabilityX());
00552 probYgood->Fill(pix->probabilityY());
00553 }
00554 }
00555 } else {
00556 LogVerbatim("TestTrackHits") << "bad hit" ;
00557 title.str("");
00558 title << "Chi2BadHit_" << subdetId << "-" << layerId;
00559 hChi2BadHit[title.str()]->Fill( chi2increment );
00560 hTotChi2BadHit->Fill( chi2increment );
00561 goodbadmerged->Fill(2);
00562 if (pix) {
00563 probXbad->Fill(pix->probabilityX());
00564 probYbad->Fill(pix->probabilityY());
00565 }
00566 }
00567 hGoodHit_vs_Chi2->Fill(chi2increment,goodhit);
00568
00569 LocalVector shitLMom;
00570 LocalPoint shitLPos;
00571 if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())) {
00572 if (simhitvecsize>2 && goodhit) {
00573 if (ioniOnly&&!shared) goodbadmerged->Fill(3);
00574 else if(!ioniOnly&&!shared) goodbadmerged->Fill(4);
00575 else goodbadmerged->Fill(5);
00576 }
00577 else if (goodhit) goodbadmerged->Fill(1);
00578 double rechitmatchedx = rhit->localPosition().x();
00579 double rechitmatchedy = rhit->localPosition().y();
00580 double mindist = 999999;
00581 float distx, disty;
00582 std::pair<LocalPoint,LocalVector> closestPair;
00583 const StripGeomDetUnit* stripDet =(StripGeomDetUnit*) ((const GluedGeomDet *)(rhit)->det())->stereoDet();
00584 const BoundPlane& plane = (rhit)->det()->surface();
00585 for(std::vector<PSimHit>::const_iterator m=assSimHits.begin(); m<assSimHits.end(); m++){
00586
00587 std::pair<LocalPoint,LocalVector> hitPair = projectHit((*m),stripDet,plane);
00588 distx = fabs(rechitmatchedx - hitPair.first.x());
00589 disty = fabs(rechitmatchedy - hitPair.first.y());
00590 double dist = distx*distx+disty*disty;
00591 if(sqrt(dist)<mindist){
00592 mindist = dist;
00593 closestPair = hitPair;
00594 }
00595 }
00596 shitLPos = closestPair.first;
00597 shitLMom = closestPair.second;
00598 } else {
00599 if (simhitvecsize>1 && goodhit) {
00600 if (ioniOnly&&!shared) goodbadmerged->Fill(3);
00601 else if(!ioniOnly&&!shared) goodbadmerged->Fill(4);
00602 else goodbadmerged->Fill(5);
00603 }
00604 else if (goodhit) goodbadmerged->Fill(1);
00605 shitLPos = shit.localPosition();
00606 shitLMom = shit.momentumAtEntry();
00607 }
00608 GlobalVector shitGMom = surf->toGlobal(shitLMom);
00609 GlobalPoint shitGPos = surf->toGlobal(shitLPos);
00610
00611 GlobalVector tsosGMom = state.globalMomentum();
00612 GlobalError tsosGMEr(state.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00613 GlobalPoint tsosGPos = state.globalPosition();
00614 GlobalError tsosGPEr = state.cartesianError().position();
00615
00616 GlobalPoint rhitGPos = (rhit)->globalPosition();
00617 GlobalError rhitGPEr = (rhit)->globalPositionError();
00618
00619 LogVerbatim("TestTrackHits") << "assSimHits.size()=" << assSimHits.size() ;
00620 LogVerbatim("TestTrackHits") << "tsos globalPos =" << tsosGPos ;
00621 LogVerbatim("TestTrackHits") << "sim hit globalPos=" << shitGPos ;
00622 LogVerbatim("TestTrackHits") << "rec hit globalPos=" << rhitGPos ;
00623 LogVerbatim("TestTrackHits") << "geographicalId =" << rhit->det()->geographicalId().rawId() ;
00624 LogVerbatim("TestTrackHits") << "surface position =" << surf->position() ;
00625
00626 # if 0
00627 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->geographicalId()="
00628 << rhit->detUnit()->geographicalId().rawId() ;
00629 LogTrace("TestTrackHits") << "rhit->det()->surface().position()="
00630 << rhit->det()->surface().position() ;
00631 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->surface().position()="
00632 << rhit->detUnit()->surface().position() ;
00633 LogTrace("TestTrackHits") << "rhit->det()->position()=" << rhit->det()->position() ;
00634 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->position()="
00635 << rhit->detUnit()->position() ;
00636 LogTrace("TestTrackHits") << "rhit->det()->surface().bounds().length()="
00637 << rhit->det()->surface().bounds().length() ;
00638 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->surface().bounds().length()="
00639 << rhit->detUnit()->surface().bounds().length() ;
00640 LogTrace("TestTrackHits") << "rhit->det()->surface().bounds().width()="
00641 << rhit->det()->surface().bounds().width() ;
00642 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->surface().bounds().width()="
00643 << rhit->detUnit()->surface().bounds().width() ;
00644 LogTrace("TestTrackHits") << "rhit->det()->surface().bounds().thickness()="
00645 << rhit->det()->surface().bounds().thickness() ;
00646 if (rhit->detUnit()) LogTrace("TestTrackHits") << "rhit->detUnit()->surface().bounds().thickness()="
00647 << rhit->detUnit()->surface().bounds().thickness() ;
00648 #endif
00649
00650 double pullGPX_rs = (rhitGPos.x()-shitGPos.x())/sqrt(rhitGPEr.cxx());
00651 double pullGPY_rs = (rhitGPos.y()-shitGPos.y())/sqrt(rhitGPEr.cyy());
00652 double pullGPZ_rs = (rhitGPos.z()-shitGPos.z())/sqrt(rhitGPEr.czz());
00653
00654
00655
00656
00657 LogTrace("TestTrackHits") << "rs" ;
00658
00659 title.str("");
00660 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs";
00661 hPullGP_X_rs[title.str()]->Fill( pullGPX_rs );
00662 title.str("");
00663 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs";
00664 hPullGP_Y_rs[title.str()]->Fill( pullGPY_rs );
00665 title.str("");
00666 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs";
00667 hPullGP_Z_rs[title.str()]->Fill( pullGPZ_rs );
00668
00669 double pullGPX_tr = (tsosGPos.x()-rhitGPos.x())/sqrt(tsosGPEr.cxx()+rhitGPEr.cxx());
00670 double pullGPY_tr = (tsosGPos.y()-rhitGPos.y())/sqrt(tsosGPEr.cyy()+rhitGPEr.cyy());
00671 double pullGPZ_tr = (tsosGPos.z()-rhitGPos.z())/sqrt(tsosGPEr.czz()+rhitGPEr.czz());
00672
00673
00674
00675
00676 LogTrace("TestTrackHits") << "tr" ;
00677
00678 title.str("");
00679 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr";
00680 hPullGP_X_tr[title.str()]->Fill( pullGPX_tr );
00681 title.str("");
00682 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr";
00683 hPullGP_Y_tr[title.str()]->Fill( pullGPY_tr );
00684 title.str("");
00685 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr";
00686 hPullGP_Z_tr[title.str()]->Fill( pullGPZ_tr );
00687
00688 double pullGPX_ts = (tsosGPos.x()-shitGPos.x())/sqrt(tsosGPEr.cxx());
00689 double pullGPY_ts = (tsosGPos.y()-shitGPos.y())/sqrt(tsosGPEr.cyy());
00690 double pullGPZ_ts = (tsosGPos.z()-shitGPos.z())/sqrt(tsosGPEr.czz());
00691
00692
00693
00694
00695 LogTrace("TestTrackHits") << "ts1" ;
00696
00697 title.str("");
00698 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts";
00699 hPullGP_X_ts[title.str()]->Fill( pullGPX_ts );
00700 title.str("");
00701 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts";
00702 hPullGP_Y_ts[title.str()]->Fill( pullGPY_ts );
00703 title.str("");
00704 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts";
00705 hPullGP_Z_ts[title.str()]->Fill( pullGPZ_ts );
00706
00707 double pullGMX_ts = (tsosGMom.x()-shitGMom.x())/sqrt(tsosGMEr.cxx());
00708 double pullGMY_ts = (tsosGMom.y()-shitGMom.y())/sqrt(tsosGMEr.cyy());
00709 double pullGMZ_ts = (tsosGMom.z()-shitGMom.z())/sqrt(tsosGMEr.czz());
00710
00711
00712
00713
00714 LogTrace("TestTrackHits") << "ts2" ;
00715
00716 title.str("");
00717 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts";
00718 hPullGM_X_ts[title.str()]->Fill( pullGMX_ts );
00719 title.str("");
00720 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts";
00721 hPullGM_Y_ts[title.str()]->Fill( pullGMY_ts );
00722 title.str("");
00723 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts";
00724 hPullGM_Z_ts[title.str()]->Fill( pullGMZ_ts );
00725 if (dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())) {
00726 Propagator* thePropagatorAnyDir = new PropagatorWithMaterial(anyDirection,0.105,theMF.product(),1.6);
00727
00728 LogTrace("TestTrackHits") << "MONO HIT" ;
00729 auto m = dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())->monoHit();
00730 CTTRHp tMonoHit = theBuilder->build(&m);
00731 if (tMonoHit==0) continue;
00732 vector<PSimHit> assMonoSimHits = hitAssociator->associateHit(*tMonoHit->hit());
00733 if (assMonoSimHits.size()==0) continue;
00734 const PSimHit sMonoHit = *(assMonoSimHits.begin());
00735 const Surface * monoSurf = &( tMonoHit->det()->surface() );
00736 if (monoSurf==0) continue;
00737 TSOS monoState = thePropagatorAnyDir->propagate(state,*monoSurf);
00738 if (monoState.isValid()==0) continue;
00739
00740 LocalVector monoShitLMom = sMonoHit.momentumAtEntry();
00741 GlobalVector monoShitGMom = monoSurf->toGlobal(monoShitLMom);
00742 LocalPoint monoShitLPos = sMonoHit.localPosition();
00743 GlobalPoint monoShitGPos = monoSurf->toGlobal(monoShitLPos);
00744
00745
00746
00747
00748 GlobalVector monoTsosGMom = monoState.globalMomentum();
00749 GlobalError monoTsosGMEr(monoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00750 GlobalPoint monoTsosGPos = monoState.globalPosition();
00751 GlobalError monoTsosGPEr = monoState.cartesianError().position();
00752
00753 GlobalPoint monoRhitGPos = tMonoHit->globalPosition();
00754 GlobalError monoRhitGPEr = tMonoHit->globalPositionError();
00755
00756 double pullGPX_rs_mono = (monoRhitGPos.x()-monoShitGPos.x())/sqrt(monoRhitGPEr.cxx());
00757 double pullGPY_rs_mono = (monoRhitGPos.y()-monoShitGPos.y())/sqrt(monoRhitGPEr.cyy());
00758 double pullGPZ_rs_mono = (monoRhitGPos.z()-monoShitGPos.z())/sqrt(monoRhitGPEr.czz());
00759
00760
00761
00762
00763 MeasurementExtractor meMo(monoState);
00764 double chi2mono = computeChi2Increment(meMo,tMonoHit);
00765
00766 title.str("");
00767 title << "Chi2Increment_mono_" << subdetId << "-" << layerId ;
00768 hChi2Increment_mono[title.str()]->Fill(chi2mono);
00769
00770 title.str("");
00771 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_mono";
00772 hPullGP_X_rs_mono[title.str()]->Fill( pullGPX_rs_mono );
00773 title.str("");
00774 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_mono";
00775 hPullGP_Y_rs_mono[title.str()]->Fill( pullGPY_rs_mono );
00776 title.str("");
00777 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_mono";
00778 hPullGP_Z_rs_mono[title.str()]->Fill( pullGPZ_rs_mono );
00779
00780 double pullGPX_tr_mono = (monoTsosGPos.x()-monoRhitGPos.x())/sqrt(monoTsosGPEr.cxx()+monoRhitGPEr.cxx());
00781 double pullGPY_tr_mono = (monoTsosGPos.y()-monoRhitGPos.y())/sqrt(monoTsosGPEr.cyy()+monoRhitGPEr.cyy());
00782 double pullGPZ_tr_mono = (monoTsosGPos.z()-monoRhitGPos.z())/sqrt(monoTsosGPEr.czz()+monoRhitGPEr.czz());
00783
00784
00785
00786
00787 title.str("");
00788 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_mono";
00789 hPullGP_X_tr_mono[title.str()]->Fill( pullGPX_tr_mono );
00790 title.str("");
00791 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_mono";
00792 hPullGP_Y_tr_mono[title.str()]->Fill( pullGPY_tr_mono );
00793 title.str("");
00794 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_mono";
00795 hPullGP_Z_tr_mono[title.str()]->Fill( pullGPZ_tr_mono );
00796
00797 double pullGPX_ts_mono = (monoTsosGPos.x()-monoShitGPos.x())/sqrt(monoTsosGPEr.cxx());
00798 double pullGPY_ts_mono = (monoTsosGPos.y()-monoShitGPos.y())/sqrt(monoTsosGPEr.cyy());
00799 double pullGPZ_ts_mono = (monoTsosGPos.z()-monoShitGPos.z())/sqrt(monoTsosGPEr.czz());
00800
00801
00802
00803
00804 title.str("");
00805 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_mono";
00806 hPullGP_X_ts_mono[title.str()]->Fill( pullGPX_ts_mono );
00807 title.str("");
00808 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_mono";
00809 hPullGP_Y_ts_mono[title.str()]->Fill( pullGPY_ts_mono );
00810 title.str("");
00811 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_mono";
00812 hPullGP_Z_ts_mono[title.str()]->Fill( pullGPZ_ts_mono );
00813
00814 double pullGMX_ts_mono = (monoTsosGMom.x()-monoShitGMom.x())/sqrt(monoTsosGMEr.cxx());
00815 double pullGMY_ts_mono = (monoTsosGMom.y()-monoShitGMom.y())/sqrt(monoTsosGMEr.cyy());
00816 double pullGMZ_ts_mono = (monoTsosGMom.z()-monoShitGMom.z())/sqrt(monoTsosGMEr.czz());
00817
00818
00819
00820
00821 title.str("");
00822 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_mono";
00823 hPullGM_X_ts_mono[title.str()]->Fill( pullGMX_ts_mono );
00824 title.str("");
00825 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_mono";
00826 hPullGM_Y_ts_mono[title.str()]->Fill( pullGMY_ts_mono );
00827 title.str("");
00828 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_mono";
00829 hPullGM_Z_ts_mono[title.str()]->Fill( pullGMZ_ts_mono );
00830
00831
00832 LogTrace("TestTrackHits") << "STEREO HIT" ;
00833 auto s = dynamic_cast<const SiStripMatchedRecHit2D*>(rhit->hit())->stereoHit();
00834 CTTRHp tStereoHit = theBuilder->build(&s);
00835 if (tStereoHit==0) continue;
00836 vector<PSimHit> assStereoSimHits = hitAssociator->associateHit(*tStereoHit->hit());
00837 if (assStereoSimHits.size()==0) continue;
00838 const PSimHit sStereoHit = *(assStereoSimHits.begin());
00839 const Surface * stereoSurf = &( tStereoHit->det()->surface() );
00840 if (stereoSurf==0) continue;
00841 TSOS stereoState = thePropagatorAnyDir->propagate(state,*stereoSurf);
00842 if (stereoState.isValid()==0) continue;
00843
00844 LocalVector stereoShitLMom = sStereoHit.momentumAtEntry();
00845 GlobalVector stereoShitGMom = stereoSurf->toGlobal(stereoShitLMom);
00846 LocalPoint stereoShitLPos = sStereoHit.localPosition();
00847 GlobalPoint stereoShitGPos = stereoSurf->toGlobal(stereoShitLPos);
00848
00849
00850
00851
00852 GlobalVector stereoTsosGMom = stereoState.globalMomentum();
00853 GlobalError stereoTsosGMEr(stereoState.cartesianError().matrix().Sub<AlgebraicSymMatrix33>(3,3));
00854 GlobalPoint stereoTsosGPos = stereoState.globalPosition();
00855 GlobalError stereoTsosGPEr = stereoState.cartesianError().position();
00856
00857 GlobalPoint stereoRhitGPos = tStereoHit->globalPosition();
00858 GlobalError stereoRhitGPEr = tStereoHit->globalPositionError();
00859
00860 MeasurementExtractor meSt(stereoState);
00861 double chi2stereo = computeChi2Increment(meSt,tStereoHit);
00862
00863 title.str("");
00864 title << "Chi2Increment_stereo_" << subdetId << "-" << layerId ;
00865 hChi2Increment_stereo[title.str()]->Fill(chi2stereo);
00866
00867 double pullGPX_rs_stereo = (stereoRhitGPos.x()-stereoShitGPos.x())/sqrt(stereoRhitGPEr.cxx());
00868 double pullGPY_rs_stereo = (stereoRhitGPos.y()-stereoShitGPos.y())/sqrt(stereoRhitGPEr.cyy());
00869 double pullGPZ_rs_stereo = (stereoRhitGPos.z()-stereoShitGPos.z())/sqrt(stereoRhitGPEr.czz());
00870
00871
00872
00873
00874 title.str("");
00875 title << "PullGP_X_" << subdetId << "-" << layerId << "_rs_stereo";
00876 hPullGP_X_rs_stereo[title.str()]->Fill( pullGPX_rs_stereo );
00877 title.str("");
00878 title << "PullGP_Y_" << subdetId << "-" << layerId << "_rs_stereo";
00879 hPullGP_Y_rs_stereo[title.str()]->Fill( pullGPY_rs_stereo );
00880 title.str("");
00881 title << "PullGP_Z_" << subdetId << "-" << layerId << "_rs_stereo";
00882 hPullGP_Z_rs_stereo[title.str()]->Fill( pullGPZ_rs_stereo );
00883
00884 double pullGPX_tr_stereo = (stereoTsosGPos.x()-stereoRhitGPos.x())/sqrt(stereoTsosGPEr.cxx()+stereoRhitGPEr.cxx());
00885 double pullGPY_tr_stereo = (stereoTsosGPos.y()-stereoRhitGPos.y())/sqrt(stereoTsosGPEr.cyy()+stereoRhitGPEr.cyy());
00886 double pullGPZ_tr_stereo = (stereoTsosGPos.z()-stereoRhitGPos.z())/sqrt(stereoTsosGPEr.czz()+stereoRhitGPEr.czz());
00887
00888
00889
00890
00891 title.str("");
00892 title << "PullGP_X_" << subdetId << "-" << layerId << "_tr_stereo";
00893 hPullGP_X_tr_stereo[title.str()]->Fill( pullGPX_tr_stereo );
00894 title.str("");
00895 title << "PullGP_Y_" << subdetId << "-" << layerId << "_tr_stereo";
00896 hPullGP_Y_tr_stereo[title.str()]->Fill( pullGPY_tr_stereo );
00897 title.str("");
00898 title << "PullGP_Z_" << subdetId << "-" << layerId << "_tr_stereo";
00899 hPullGP_Z_tr_stereo[title.str()]->Fill( pullGPZ_tr_stereo );
00900
00901 double pullGPX_ts_stereo = (stereoTsosGPos.x()-stereoShitGPos.x())/sqrt(stereoTsosGPEr.cxx());
00902 double pullGPY_ts_stereo = (stereoTsosGPos.y()-stereoShitGPos.y())/sqrt(stereoTsosGPEr.cyy());
00903 double pullGPZ_ts_stereo = (stereoTsosGPos.z()-stereoShitGPos.z())/sqrt(stereoTsosGPEr.czz());
00904
00905
00906
00907
00908 title.str("");
00909 title << "PullGP_X_" << subdetId << "-" << layerId << "_ts_stereo";
00910 hPullGP_X_ts_stereo[title.str()]->Fill( pullGPX_ts_stereo );
00911 title.str("");
00912 title << "PullGP_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00913 hPullGP_Y_ts_stereo[title.str()]->Fill( pullGPY_ts_stereo );
00914 title.str("");
00915 title << "PullGP_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00916 hPullGP_Z_ts_stereo[title.str()]->Fill( pullGPZ_ts_stereo );
00917
00918 double pullGMX_ts_stereo = (stereoTsosGMom.x()-stereoShitGMom.x())/sqrt(stereoTsosGMEr.cxx());
00919 double pullGMY_ts_stereo = (stereoTsosGMom.y()-stereoShitGMom.y())/sqrt(stereoTsosGMEr.cyy());
00920 double pullGMZ_ts_stereo = (stereoTsosGMom.z()-stereoShitGMom.z())/sqrt(stereoTsosGMEr.czz());
00921
00922
00923
00924
00925 title.str("");
00926 title << "PullGM_X_" << subdetId << "-" << layerId << "_ts_stereo";
00927 hPullGM_X_ts_stereo[title.str()]->Fill( pullGMX_ts_stereo );
00928 title.str("");
00929 title << "PullGM_Y_" << subdetId << "-" << layerId << "_ts_stereo";
00930 hPullGM_Y_ts_stereo[title.str()]->Fill( pullGMY_ts_stereo );
00931 title.str("");
00932 title << "PullGM_Z_" << subdetId << "-" << layerId << "_ts_stereo";
00933 hPullGM_Z_ts_stereo[title.str()]->Fill( pullGMZ_ts_stereo );
00934 }
00935 }
00936 LogTrace("TestTrackHits") << "traj chi2=" << tchi2 ;
00937 LogTrace("TestTrackHits") << "track chi2=" << track->chi2() ;
00938 i++;
00939 }
00940 LogTrace("TestTrackHits") << "end of event: processd hits=" << evtHits ;
00941 delete hitAssociator;
00942 }
00943
00944 void TestTrackHits::endJob() {
00945
00946 TDirectory * chi2d = file->mkdir("Chi2Increment");
00947
00948 TDirectory * gp_ts = file->mkdir("GP_TSOS-SimHit");
00949 TDirectory * gm_ts = file->mkdir("GM_TSOS-SimHit");
00950 TDirectory * gp_tr = file->mkdir("GP_TSOS-RecHit");
00951 TDirectory * gp_rs = file->mkdir("GP_RecHit-SimHit");
00952
00953 TDirectory * gp_tsx = gp_ts->mkdir("X");
00954 TDirectory * gp_tsy = gp_ts->mkdir("Y");
00955 TDirectory * gp_tsz = gp_ts->mkdir("Z");
00956 TDirectory * gm_tsx = gm_ts->mkdir("X");
00957 TDirectory * gm_tsy = gm_ts->mkdir("Y");
00958 TDirectory * gm_tsz = gm_ts->mkdir("Z");
00959 TDirectory * gp_trx = gp_tr->mkdir("X");
00960 TDirectory * gp_try = gp_tr->mkdir("Y");
00961 TDirectory * gp_trz = gp_tr->mkdir("Z");
00962 TDirectory * gp_rsx = gp_rs->mkdir("X");
00963 TDirectory * gp_rsy = gp_rs->mkdir("Y");
00964 TDirectory * gp_rsz = gp_rs->mkdir("Z");
00965
00966 TDirectory * gp_tsx_mono = gp_ts->mkdir("MONOX");
00967 TDirectory * gp_tsy_mono = gp_ts->mkdir("MONOY");
00968 TDirectory * gp_tsz_mono = gp_ts->mkdir("MONOZ");
00969 TDirectory * gm_tsx_mono = gm_ts->mkdir("MONOX");
00970 TDirectory * gm_tsy_mono = gm_ts->mkdir("MONOY");
00971 TDirectory * gm_tsz_mono = gm_ts->mkdir("MONOZ");
00972 TDirectory * gp_trx_mono = gp_tr->mkdir("MONOX");
00973 TDirectory * gp_try_mono = gp_tr->mkdir("MONOY");
00974 TDirectory * gp_trz_mono = gp_tr->mkdir("MONOZ");
00975 TDirectory * gp_rsx_mono = gp_rs->mkdir("MONOX");
00976 TDirectory * gp_rsy_mono = gp_rs->mkdir("MONOY");
00977 TDirectory * gp_rsz_mono = gp_rs->mkdir("MONOZ");
00978
00979 TDirectory * gp_tsx_stereo = gp_ts->mkdir("STEREOX");
00980 TDirectory * gp_tsy_stereo = gp_ts->mkdir("STEREOY");
00981 TDirectory * gp_tsz_stereo = gp_ts->mkdir("STEREOZ");
00982 TDirectory * gm_tsx_stereo = gm_ts->mkdir("STEREOX");
00983 TDirectory * gm_tsy_stereo = gm_ts->mkdir("STEREOY");
00984 TDirectory * gm_tsz_stereo = gm_ts->mkdir("STEREOZ");
00985 TDirectory * gp_trx_stereo = gp_tr->mkdir("STEREOX");
00986 TDirectory * gp_try_stereo = gp_tr->mkdir("STEREOY");
00987 TDirectory * gp_trz_stereo = gp_tr->mkdir("STEREOZ");
00988 TDirectory * gp_rsx_stereo = gp_rs->mkdir("STEREOX");
00989 TDirectory * gp_rsy_stereo = gp_rs->mkdir("STEREOY");
00990 TDirectory * gp_rsz_stereo = gp_rs->mkdir("STEREOZ");
00991
00992 chi2d->cd();
00993 hTotChi2Increment->Write();
00994 hTotChi2GoodHit->Write();
00995 hTotChi2BadHit->Write();
00996 hTotChi2DeltaHit->Write();
00997 hTotChi2NSharedHit->Write();
00998 hTotChi2SharedHit->Write();
00999 hProcess_vs_Chi2->Write();
01000 hClsize_vs_Chi2->Write();
01001 hPixClsize_vs_Chi2->Write();
01002 hPrjClsize_vs_Chi2->Write();
01003 hSt1Clsize_vs_Chi2->Write();
01004 hSt2Clsize_vs_Chi2->Write();
01005 hGoodHit_vs_Chi2->Write();
01006 hClusterSize->Write();
01007 hPixClusterSize->Write();
01008 hPrjClusterSize->Write();
01009 hSt1ClusterSize->Write();
01010 hSt2ClusterSize->Write();
01011 hSimHitVecSize->Write();
01012 hPixSimHitVecSize->Write();
01013 hPrjSimHitVecSize->Write();
01014 hSt1SimHitVecSize->Write();
01015 hSt2SimHitVecSize->Write();
01016 goodbadmerged->Write();
01017 energyLossRatio->Write();
01018 mergedPull->Write();
01019 probXgood->Write();
01020 probXbad->Write();
01021 probXdelta->Write();
01022 probXshared->Write();
01023 probXnoshare->Write();
01024 probYgood->Write();
01025 probYbad->Write();
01026 probYdelta->Write();
01027 probYshared->Write();
01028 probYnoshare->Write();
01029 for (int i=0; i!=6; i++)
01030 for (int j=0; j!=9; j++){
01031 if (i==0 && j>2) break;
01032 if (i==1 && j>1) break;
01033 if (i==2 && j>3) break;
01034 if (i==3 && j>2) break;
01035 if (i==4 && j>5) break;
01036 if (i==5 && j>8) break;
01037 chi2d->cd();
01038 title.str("");
01039 title << "Chi2Increment_" << i+1 << "-" << j+1 ;
01040 hChi2Increment[title.str()]->Write();
01041 title.str("");
01042 title << "Chi2IncrementVsEta_" << i+1 << "-" << j+1 ;
01043 hChi2IncrementVsEta[title.str()]->Write();
01044 title.str("");
01045 title << "Chi2GoodHit_" << i+1 << "-" << j+1 ;
01046 hChi2GoodHit[title.str()]->Write();
01047 title.str("");
01048 title << "Chi2BadHit_" << i+1 << "-" << j+1 ;
01049 hChi2BadHit[title.str()]->Write();
01050 title.str("");
01051 title << "Chi2DeltaHit_" << i+1 << "-" << j+1 ;
01052 hChi2DeltaHit[title.str()]->Write();
01053 title.str("");
01054 title << "Chi2NSharedHit_" << i+1 << "-" << j+1 ;
01055 hChi2NSharedHit[title.str()]->Write();
01056 title.str("");
01057 title << "Chi2SharedHit_" << i+1 << "-" << j+1 ;
01058 hChi2SharedHit[title.str()]->Write();
01059
01060 gp_ts->cd();
01061 gp_tsx->cd();
01062 title.str("");
01063 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts";
01064 hPullGP_X_ts[title.str()]->Write();
01065 gp_tsy->cd();
01066 title.str("");
01067 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts";
01068 hPullGP_Y_ts[title.str()]->Write();
01069 gp_tsz->cd();
01070 title.str("");
01071 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts";
01072 hPullGP_Z_ts[title.str()]->Write();
01073
01074 gm_ts->cd();
01075 gm_tsx->cd();
01076 title.str("");
01077 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts";
01078 hPullGM_X_ts[title.str()]->Write();
01079 gm_tsy->cd();
01080 title.str("");
01081 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts";
01082 hPullGM_Y_ts[title.str()]->Write();
01083 gm_tsz->cd();
01084 title.str("");
01085 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts";
01086 hPullGM_Z_ts[title.str()]->Write();
01087
01088 gp_tr->cd();
01089 gp_trx->cd();
01090 title.str("");
01091 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr";
01092 hPullGP_X_tr[title.str()]->Write();
01093 gp_try->cd();
01094 title.str("");
01095 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr";
01096 hPullGP_Y_tr[title.str()]->Write();
01097 gp_trz->cd();
01098 title.str("");
01099 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr";
01100 hPullGP_Z_tr[title.str()]->Write();
01101
01102 gp_rs->cd();
01103 gp_rsx->cd();
01104 title.str("");
01105 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs";
01106 hPullGP_X_rs[title.str()]->Write();
01107 gp_rsy->cd();
01108 title.str("");
01109 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs";
01110 hPullGP_Y_rs[title.str()]->Write();
01111 gp_rsz->cd();
01112 title.str("");
01113 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs";
01114 hPullGP_Z_rs[title.str()]->Write();
01115
01116 if ( ((i==2||i==4)&&(j==0||j==1)) || (i==3||i==5) ){
01117 chi2d->cd();
01118 title.str("");
01119 title << "Chi2Increment_mono_" << i+1 << "-" << j+1 ;
01120 hChi2Increment_mono[title.str()]->Write();
01121 title.str("");
01122 title << "Chi2Increment_stereo_" << i+1 << "-" << j+1 ;
01123 hChi2Increment_stereo[title.str()]->Write();
01124
01125 gp_ts->cd();
01126 gp_tsx_mono->cd();
01127 title.str("");
01128 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_mono";
01129 hPullGP_X_ts_mono[title.str()]->Write();
01130 gp_tsy_mono->cd();
01131 title.str("");
01132 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_mono";
01133 hPullGP_Y_ts_mono[title.str()]->Write();
01134 gp_tsz_mono->cd();
01135 title.str("");
01136 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_mono";
01137 hPullGP_Z_ts_mono[title.str()]->Write();
01138
01139 gm_ts->cd();
01140 gm_tsx_mono->cd();
01141 title.str("");
01142 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_mono";
01143 hPullGM_X_ts_mono[title.str()]->Write();
01144 gm_tsy_mono->cd();
01145 title.str("");
01146 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_mono";
01147 hPullGM_Y_ts_mono[title.str()]->Write();
01148 gm_tsz_mono->cd();
01149 title.str("");
01150 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_mono";
01151 hPullGM_Z_ts_mono[title.str()]->Write();
01152
01153 gp_tr->cd();
01154 gp_trx_mono->cd();
01155 title.str("");
01156 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_mono";
01157 hPullGP_X_tr_mono[title.str()]->Write();
01158 gp_try_mono->cd();
01159 title.str("");
01160 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_mono";
01161 hPullGP_Y_tr_mono[title.str()]->Write();
01162 gp_trz_mono->cd();
01163 title.str("");
01164 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_mono";
01165 hPullGP_Z_tr_mono[title.str()]->Write();
01166
01167 gp_rs->cd();
01168 gp_rsx_mono->cd();
01169 title.str("");
01170 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_mono";
01171 hPullGP_X_rs_mono[title.str()]->Write();
01172 gp_rsy_mono->cd();
01173 title.str("");
01174 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_mono";
01175 hPullGP_Y_rs_mono[title.str()]->Write();
01176 gp_rsz_mono->cd();
01177 title.str("");
01178 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_mono";
01179 hPullGP_Z_rs_mono[title.str()]->Write();
01180
01181
01182 gp_ts->cd();
01183 gp_tsx_stereo->cd();
01184 title.str("");
01185 title << "PullGP_X_" << i+1 << "-" << j+1 << "_ts_stereo";
01186 hPullGP_X_ts_stereo[title.str()]->Write();
01187 gp_tsy_stereo->cd();
01188 title.str("");
01189 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
01190 hPullGP_Y_ts_stereo[title.str()]->Write();
01191 gp_tsz_stereo->cd();
01192 title.str("");
01193 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
01194 hPullGP_Z_ts_stereo[title.str()]->Write();
01195
01196 gm_ts->cd();
01197 gm_tsx_stereo->cd();
01198 title.str("");
01199 title << "PullGM_X_" << i+1 << "-" << j+1 << "_ts_stereo";
01200 hPullGM_X_ts_stereo[title.str()]->Write();
01201 gm_tsy_stereo->cd();
01202 title.str("");
01203 title << "PullGM_Y_" << i+1 << "-" << j+1 << "_ts_stereo";
01204 hPullGM_Y_ts_stereo[title.str()]->Write();
01205 gm_tsz_stereo->cd();
01206 title.str("");
01207 title << "PullGM_Z_" << i+1 << "-" << j+1 << "_ts_stereo";
01208 hPullGM_Z_ts_stereo[title.str()]->Write();
01209
01210 gp_tr->cd();
01211 gp_trx_stereo->cd();
01212 title.str("");
01213 title << "PullGP_X_" << i+1 << "-" << j+1 << "_tr_stereo";
01214 hPullGP_X_tr_stereo[title.str()]->Write();
01215 gp_try_stereo->cd();
01216 title.str("");
01217 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_tr_stereo";
01218 hPullGP_Y_tr_stereo[title.str()]->Write();
01219 gp_trz_stereo->cd();
01220 title.str("");
01221 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_tr_stereo";
01222 hPullGP_Z_tr_stereo[title.str()]->Write();
01223
01224 gp_rs->cd();
01225 gp_rsx_stereo->cd();
01226 title.str("");
01227 title << "PullGP_X_" << i+1 << "-" << j+1 << "_rs_stereo";
01228 hPullGP_X_rs_stereo[title.str()]->Write();
01229 gp_rsy_stereo->cd();
01230 title.str("");
01231 title << "PullGP_Y_" << i+1 << "-" << j+1 << "_rs_stereo";
01232 hPullGP_Y_rs_stereo[title.str()]->Write();
01233 gp_rsz_stereo->cd();
01234 title.str("");
01235 title << "PullGP_Z_" << i+1 << "-" << j+1 << "_rs_stereo";
01236 hPullGP_Z_rs_stereo[title.str()]->Write();
01237 }
01238 }
01239
01240 file->Close();
01241 }
01242
01243
01244
01245
01246 std::pair<LocalPoint,LocalVector>
01247 TestTrackHits::projectHit( const PSimHit& hit, const StripGeomDetUnit* stripDet, const BoundPlane& plane)
01248 {
01249 const StripTopology& topol = stripDet->specificTopology();
01250 GlobalPoint globalpos= stripDet->surface().toGlobal(hit.localPosition());
01251 LocalPoint localHit = plane.toLocal(globalpos);
01252
01253 LocalVector locdir=hit.localDirection();
01254
01255
01256 GlobalVector globaldir= stripDet->surface().toGlobal(locdir);
01257 LocalVector dir=plane.toLocal(globaldir);
01258 float scale = -localHit.z() / dir.z();
01259
01260 LocalPoint projectedPos = localHit + scale*dir;
01261
01262 float selfAngle = topol.stripAngle( topol.strip( hit.localPosition()));
01263
01264 LocalVector stripDir( sin(selfAngle), cos(selfAngle), 0);
01265
01266 LocalVector localStripDir( plane.toLocal(stripDet->surface().toGlobal( stripDir)));
01267
01268 return std::pair<LocalPoint,LocalVector>( projectedPos, localStripDir);
01269 }
01270
01271 template<unsigned int D>
01272 double TestTrackHits::computeChi2Increment(MeasurementExtractor me,
01273 TransientTrackingRecHit::ConstRecHitPointer rhit) {
01274 typedef typename AlgebraicROOTObject<D>::Vector VecD;
01275 typedef typename AlgebraicROOTObject<D,D>::SymMatrix SMatDD;
01276 VecD r = asSVector<D>(rhit->parameters()) - me.measuredParameters<D>(*rhit);
01277
01278 SMatDD R = asSMatrix<D>(rhit->parametersError()) + me.measuredError<D>(*rhit);
01279 R.Invert();
01280 return ROOT::Math::Similarity(r,R) ;
01281 }
01282
01283 #include "FWCore/Framework/interface/ModuleFactory.h"
01284 #include "FWCore/Framework/interface/MakerMacros.h"
01285 DEFINE_FWK_MODULE(TestTrackHits);