CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_10_patch2/src/DQM/TrackingMonitor/src/TrackBuildingAnalyzer.cc

Go to the documentation of this file.
00001 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateTransform.h"
00002 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
00003 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
00004 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
00005 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00006 #include "TrackingTools/PatternTools/interface/TSCBLBuilderNoMaterial.h"
00007 #include "TrackingTools/PatternTools/interface/TSCPBuilderNoMaterial.h"
00008 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00009 #include "MagneticField/Engine/interface/MagneticField.h"
00010 
00011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00012 #include "FWCore/Utilities/interface/InputTag.h"
00013 #include "DQMServices/Core/interface/DQMStore.h"
00014 #include "DQM/TrackingMonitor/interface/TrackBuildingAnalyzer.h"
00015 #include <string>
00016 #include "TMath.h"
00017 
00018 #include <iostream>
00019 
00020 TrackBuildingAnalyzer::TrackBuildingAnalyzer(const edm::ParameterSet& iConfig) 
00021     : conf_( iConfig )
00022     , SeedPt(NULL)
00023     , SeedEta(NULL)
00024     , SeedPhi(NULL)
00025     , SeedPhiVsEta(NULL)
00026     , SeedTheta(NULL)
00027     , SeedQ(NULL)
00028     , SeedDxy(NULL)
00029     , SeedDz(NULL)
00030     , NumberOfRecHitsPerSeed(NULL)
00031     , NumberOfRecHitsPerSeedVsPhiProfile(NULL)
00032     , NumberOfRecHitsPerSeedVsEtaProfile(NULL)
00033 {
00034 }
00035 
00036 TrackBuildingAnalyzer::~TrackBuildingAnalyzer() 
00037 { 
00038 }
00039 
00040 void TrackBuildingAnalyzer::beginJob(DQMStore * dqmStore_) 
00041 {
00042 
00043     // parameters from the configuration
00044     std::string AlgoName       = conf_.getParameter<std::string>("AlgoName");
00045     std::string MEFolderName   = conf_.getParameter<std::string>("FolderName"); 
00046 
00047     // use the AlgoName and Quality Name 
00048     std::string CatagoryName = AlgoName;
00049 
00050     // get binning from the configuration
00051     int    TrackPtBin = conf_.getParameter<int>(   "TrackPtBin");
00052     double TrackPtMin = conf_.getParameter<double>("TrackPtMin");
00053     double TrackPtMax = conf_.getParameter<double>("TrackPtMax");
00054 
00055     int    PhiBin     = conf_.getParameter<int>(   "PhiBin");
00056     double PhiMin     = conf_.getParameter<double>("PhiMin");
00057     double PhiMax     = conf_.getParameter<double>("PhiMax");
00058 
00059     int    EtaBin     = conf_.getParameter<int>(   "EtaBin");
00060     double EtaMin     = conf_.getParameter<double>("EtaMin");
00061     double EtaMax     = conf_.getParameter<double>("EtaMax");
00062 
00063     int    ThetaBin   = conf_.getParameter<int>(   "ThetaBin");
00064     double ThetaMin   = conf_.getParameter<double>("ThetaMin");
00065     double ThetaMax   = conf_.getParameter<double>("ThetaMax");
00066 
00067     int    TrackQBin  = conf_.getParameter<int>(   "TrackQBin");
00068     double TrackQMin  = conf_.getParameter<double>("TrackQMin");
00069     double TrackQMax  = conf_.getParameter<double>("TrackQMax");
00070 
00071     int    SeedDxyBin = conf_.getParameter<int>(   "SeedDxyBin");
00072     double SeedDxyMin = conf_.getParameter<double>("SeedDxyMin");
00073     double SeedDxyMax = conf_.getParameter<double>("SeedDxyMax");
00074 
00075     int    SeedDzBin  = conf_.getParameter<int>(   "SeedDzBin");
00076     double SeedDzMin  = conf_.getParameter<double>("SeedDzMin");
00077     double SeedDzMax  = conf_.getParameter<double>("SeedDzMax");
00078 
00079     int    SeedHitBin = conf_.getParameter<int>(   "SeedHitBin");
00080     double SeedHitMin = conf_.getParameter<double>("SeedHitMin");
00081     double SeedHitMax = conf_.getParameter<double>("SeedHitMax");
00082 
00083     int    TCDxyBin   = conf_.getParameter<int>(   "TCDxyBin");
00084     double TCDxyMin   = conf_.getParameter<double>("TCDxyMin");
00085     double TCDxyMax   = conf_.getParameter<double>("TCDxyMax");
00086 
00087     int    TCDzBin    = conf_.getParameter<int>(   "TCDzBin");
00088     double TCDzMin    = conf_.getParameter<double>("TCDzMin");
00089     double TCDzMax    = conf_.getParameter<double>("TCDzMax");
00090 
00091     int    TCHitBin   = conf_.getParameter<int>(   "TCHitBin");
00092     double TCHitMin   = conf_.getParameter<double>("TCHitMin");
00093     double TCHitMax   = conf_.getParameter<double>("TCHitMax");
00094 
00095 
00096     edm::InputTag seedProducer   = conf_.getParameter<edm::InputTag>("SeedProducer");
00097     edm::InputTag tcProducer     = conf_.getParameter<edm::InputTag>("TCProducer");
00098 
00099     doAllPlots=conf_.getParameter<bool>("doAllPlots");
00100     doAllSeedPlots=conf_.getParameter<bool>("doSeedParameterHistos");
00101     doTCPlots=conf_.getParameter<bool>("doTrackCandHistos");
00102     doPT=conf_.getParameter<bool>("doSeedPTHisto");
00103     doETA=conf_.getParameter<bool>("doSeedETAHisto");
00104     doPHI=conf_.getParameter<bool>("doSeedPHIHisto");
00105     doPHIVsETA=conf_.getParameter<bool>("doSeedPHIVsETAHisto");
00106     doTheta=conf_.getParameter<bool>("doSeedThetaHisto");
00107     doQ=conf_.getParameter<bool>("doSeedQHisto");
00108     doDxy =conf_.getParameter<bool>("doSeedDxyHisto");
00109     doDz =conf_.getParameter<bool>("doSeedDzHisto");
00110     doNRecHits=conf_.getParameter<bool>("doSeedNRecHitsHisto");
00111     doProfPHI=conf_.getParameter<bool>("doSeedNVsPhiProf");
00112     doProfETA=conf_.getParameter<bool>("doSeedNVsEtaProf");
00113 
00114     //    if (doAllPlots){doAllSeedPlots=true; doTCPlots=true;}
00115 
00116     dqmStore_->setCurrentFolder(MEFolderName);
00117 
00118     // book the Seed histograms
00119     // ---------------------------------------------------------------------------------//
00120     dqmStore_->setCurrentFolder(MEFolderName+"/TrackBuilding");
00121     
00122     if (doAllSeedPlots || doPT) {
00123       histname = "SeedPt_"+seedProducer.label() + "_";
00124       SeedPt = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TrackPtBin, TrackPtMin, TrackPtMax);
00125       SeedPt->setAxisTitle("Seed p_{T} (GeV/c)", 1);
00126       SeedPt->setAxisTitle("Number of Seeds", 2);
00127     }
00128 
00129     if (doAllSeedPlots || doETA) {
00130       histname = "SeedEta_"+seedProducer.label() + "_";
00131       SeedEta = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax);
00132       SeedEta->setAxisTitle("Seed #eta", 1);
00133       SeedEta->setAxisTitle("Number of Seeds", 2);
00134     }
00135 
00136     if (doAllSeedPlots || doPHI) {
00137       histname = "SeedPhi_"+seedProducer.label() + "_";
00138       SeedPhi = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax);
00139       SeedPhi->setAxisTitle("Seed #phi", 1);
00140       SeedPhi->setAxisTitle("Number of Seed", 2);
00141     }
00142 
00143     if (doAllSeedPlots || doPHIVsETA) {
00144       histname = "SeedPhiVsEta_"+seedProducer.label() + "_";
00145       SeedPhiVsEta = dqmStore_->book2D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, PhiBin, PhiMin, PhiMax);
00146       SeedPhiVsEta->setAxisTitle("Seed #eta", 1);
00147       SeedPhiVsEta->setAxisTitle("Seed #phi", 2);
00148     }
00149 
00150     if (doAllSeedPlots || doTheta){
00151       histname = "SeedTheta_"+seedProducer.label() + "_";
00152       SeedTheta = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, ThetaBin, ThetaMin, ThetaMax);
00153       SeedTheta->setAxisTitle("Seed #theta", 1);
00154       SeedTheta->setAxisTitle("Number of Seeds", 2);
00155     }
00156 
00157     if (doAllSeedPlots || doQ){
00158     histname = "SeedQ_"+seedProducer.label() + "_";
00159     SeedQ = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TrackQBin, TrackQMin, TrackQMax);
00160     SeedQ->setAxisTitle("Seed Charge", 1);
00161     SeedQ->setAxisTitle("Number of Seeds",2);
00162     }
00163 
00164     if (doAllSeedPlots || doDxy){
00165       histname = "SeedDxy_"+seedProducer.label() + "_";
00166       SeedDxy = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, SeedDxyBin, SeedDxyMin, SeedDxyMax);
00167       SeedDxy->setAxisTitle("Seed d_{xy} (cm)", 1);
00168       SeedDxy->setAxisTitle("Number of Seeds",2);
00169     }
00170 
00171     if (doAllSeedPlots || doDz){
00172       histname = "SeedDz_"+seedProducer.label() + "_";
00173       SeedDz = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, SeedDzBin, SeedDzMin, SeedDzMax);
00174       SeedDz->setAxisTitle("Seed d_{z} (cm)", 1);
00175       SeedDz->setAxisTitle("Number of Seeds",2);
00176     }
00177 
00178     if (doAllSeedPlots || doNRecHits){
00179       histname = "NumberOfRecHitsPerSeed_"+seedProducer.label() + "_";
00180       NumberOfRecHitsPerSeed = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, SeedHitBin, SeedHitMin, SeedHitMax);
00181       NumberOfRecHitsPerSeed->setAxisTitle("Number of RecHits per Seed", 1);
00182       NumberOfRecHitsPerSeed->setAxisTitle("Number of Seeds",2);
00183     }
00184 
00185     if (doAllSeedPlots || doProfPHI){
00186       histname = "NumberOfRecHitsPerSeedVsPhiProfile_"+seedProducer.label() + "_";
00187       NumberOfRecHitsPerSeedVsPhiProfile = dqmStore_->bookProfile(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax, SeedHitBin, SeedHitMin, SeedHitMax,"s");
00188       NumberOfRecHitsPerSeedVsPhiProfile->setAxisTitle("Seed #phi",1);
00189       NumberOfRecHitsPerSeedVsPhiProfile->setAxisTitle("Number of RecHits of each Seed",2);
00190     }
00191 
00192     if (doAllSeedPlots || doProfETA){
00193       histname = "NumberOfRecHitsPerSeedVsEtaProfile_"+seedProducer.label() + "_";
00194       NumberOfRecHitsPerSeedVsEtaProfile = dqmStore_->bookProfile(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, SeedHitBin, SeedHitMin, SeedHitMax,"s");
00195       NumberOfRecHitsPerSeedVsEtaProfile->setAxisTitle("Seed #eta",1);
00196       NumberOfRecHitsPerSeedVsEtaProfile->setAxisTitle("Number of RecHits of each Seed",2);
00197     }
00198 
00199     // book the TrackCandidate histograms
00200     // ---------------------------------------------------------------------------------//
00201 
00202     if (doTCPlots){
00203 
00204       dqmStore_->setCurrentFolder(MEFolderName+"/TrackBuilding");
00205 
00206       histname = "TrackCandPt_"+tcProducer.label() + "_";
00207       TrackCandPt = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TrackPtBin, TrackPtMin, TrackPtMax);
00208       TrackCandPt->setAxisTitle("Track Candidate p_{T} (GeV/c)", 1);
00209       TrackCandPt->setAxisTitle("Number of Track Candidates", 2);
00210       
00211       histname = "TrackCandEta_"+tcProducer.label() + "_";
00212       TrackCandEta = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax);
00213       TrackCandEta->setAxisTitle("Track Candidate #eta", 1);
00214       TrackCandEta->setAxisTitle("Number of Track Candidates", 2);
00215       
00216       histname = "TrackCandPhi_"+tcProducer.label() + "_";
00217       TrackCandPhi = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax);
00218       TrackCandPhi->setAxisTitle("Track Candidate #phi", 1);
00219       TrackCandPhi->setAxisTitle("Number of Track Candidates", 2);
00220       
00221       histname = "TrackCandTheta_"+tcProducer.label() + "_";
00222       TrackCandTheta = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, ThetaBin, ThetaMin, ThetaMax);
00223       TrackCandTheta->setAxisTitle("Track Candidate #theta", 1);
00224       TrackCandTheta->setAxisTitle("Number of Track Candidates", 2);
00225       
00226       histname = "TrackCandQ_"+tcProducer.label() + "_";
00227       TrackCandQ = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TrackQBin, TrackQMin, TrackQMax);
00228       TrackCandQ->setAxisTitle("Track Candidate Charge", 1);
00229       TrackCandQ->setAxisTitle("Number of Track Candidates",2);
00230       
00231       histname = "TrackCandDxy_"+tcProducer.label() + "_";
00232       TrackCandDxy = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TCDxyBin, TCDxyMin, TCDxyMax);
00233       TrackCandDxy->setAxisTitle("Track Candidate d_{xy} (cm)", 1);
00234       TrackCandDxy->setAxisTitle("Number of Track Candidates",2);
00235       
00236       histname = "TrackCandDz_"+tcProducer.label() + "_";
00237       TrackCandDz = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TCDzBin, TCDzMin, TCDzMax);
00238       TrackCandDz->setAxisTitle("Track Candidate d_{z} (cm)", 1);
00239       TrackCandDz->setAxisTitle("Number of Track Candidates",2);
00240       
00241       histname = "NumberOfRecHitsPerTrackCand_"+tcProducer.label() + "_";
00242       NumberOfRecHitsPerTrackCand = dqmStore_->book1D(histname+CatagoryName, histname+CatagoryName, TCHitBin, TCHitMin, TCHitMax);
00243       NumberOfRecHitsPerTrackCand->setAxisTitle("Number of RecHits per Track Candidate", 1);
00244       NumberOfRecHitsPerTrackCand->setAxisTitle("Number of Track Candidates",2);
00245       
00246       histname = "NumberOfRecHitsPerTrackCandVsPhiProfile_"+tcProducer.label() + "_";
00247       NumberOfRecHitsPerTrackCandVsPhiProfile = dqmStore_->bookProfile(histname+CatagoryName, histname+CatagoryName, PhiBin, PhiMin, PhiMax, TCHitBin, TCHitMin, TCHitMax,"s");
00248       NumberOfRecHitsPerTrackCandVsPhiProfile->setAxisTitle("Track Candidate #phi",1);
00249       NumberOfRecHitsPerTrackCandVsPhiProfile->setAxisTitle("Number of RecHits of each Track Candidate",2);
00250       
00251       histname = "NumberOfRecHitsPerTrackCandVsEtaProfile_"+tcProducer.label() + "_";
00252       NumberOfRecHitsPerTrackCandVsEtaProfile = dqmStore_->bookProfile(histname+CatagoryName, histname+CatagoryName, EtaBin, EtaMin, EtaMax, TCHitBin, TCHitMin, TCHitMax,"s");
00253       NumberOfRecHitsPerTrackCandVsEtaProfile->setAxisTitle("Track Candidate #eta",1);
00254       NumberOfRecHitsPerTrackCandVsEtaProfile->setAxisTitle("Number of RecHits of each Track Candidate",2);
00255     }
00256 }
00257 
00258 // -- Analyse
00259 // ---------------------------------------------------------------------------------//
00260 void TrackBuildingAnalyzer::analyze
00261 (
00262     const edm::Event& iEvent,
00263     const edm::EventSetup& iSetup,
00264     const TrajectorySeed& candidate,
00265     const reco::BeamSpot& bs,
00266     const edm::ESHandle<MagneticField>& theMF,
00267     const edm::ESHandle<TransientTrackingRecHitBuilder>& theTTRHBuilder
00268 )
00269 {
00270 
00271     
00272     TSCBLBuilderNoMaterial tscblBuilder;
00273 
00274     //get parameters and errors from the candidate state
00275     TransientTrackingRecHit::RecHitPointer recHit = theTTRHBuilder->build(&*(candidate.recHits().second-1));
00276     TrajectoryStateOnSurface state = trajectoryStateTransform::transientState( candidate.startingState(), recHit->surface(), theMF.product());
00277     TrajectoryStateClosestToBeamLine tsAtClosestApproachSeed = tscblBuilder(*state.freeState(),bs);//as in TrackProducerAlgorithm
00278     if(!(tsAtClosestApproachSeed.isValid())) {
00279         edm::LogVerbatim("TrackBuilding") << "TrajectoryStateClosestToBeamLine not valid";
00280         return;
00281     }
00282     GlobalPoint  v0 = tsAtClosestApproachSeed.trackStateAtPCA().position();
00283     GlobalVector p = tsAtClosestApproachSeed.trackStateAtPCA().momentum();
00284     GlobalPoint  v(v0.x()-bs.x0(),v0.y()-bs.y0(),v0.z()-bs.z0());
00285 
00286     double pt           = sqrt(state.globalMomentum().perp2());
00287     double eta          = state.globalMomentum().eta();
00288     double phi          = state.globalMomentum().phi();
00289     double theta        = state.globalMomentum().theta();
00290     //double pm           = sqrt(state.globalMomentum().mag2());
00291     //double pz           = state.globalMomentum().z();
00292     //double qoverp       = tsAtClosestApproachSeed.trackStateAtPCA().charge()/p.mag();
00293     //double theta        = p.theta();
00294     //double lambda       = M_PI/2-p.theta();
00295     double numberOfHits = candidate.recHits().second-candidate.recHits().first;
00296     double dxy          = (-v.x()*sin(p.phi())+v.y()*cos(p.phi()));
00297     double dz           = v.z() - (v.x()*p.x()+v.y()*p.y())/p.perp() * p.z()/p.perp();
00298 
00299     // fill the ME's
00300     if (doAllSeedPlots || doQ)SeedQ->Fill( state.charge() );
00301     if (doAllSeedPlots || doPT) SeedPt->Fill( pt );
00302     if (doAllSeedPlots || doETA) SeedEta->Fill( eta );
00303     if (doAllSeedPlots || doPHI) SeedPhi->Fill( phi );
00304     if (doAllSeedPlots || doPHIVsETA) SeedPhiVsEta->Fill( eta, phi);
00305     if (doAllSeedPlots || doTheta) SeedTheta->Fill( theta );
00306     if (doAllSeedPlots || doDxy) SeedDxy->Fill( dxy );
00307     if (doAllSeedPlots || doDz) SeedDz->Fill( dz );
00308     if (doAllSeedPlots || doNRecHits) NumberOfRecHitsPerSeed->Fill( numberOfHits );
00309     if (doAllSeedPlots || doProfETA) NumberOfRecHitsPerSeedVsEtaProfile->Fill( eta, numberOfHits );
00310     if (doAllSeedPlots || doProfPHI) NumberOfRecHitsPerSeedVsPhiProfile->Fill( phi, numberOfHits );
00311 }
00312 
00313 // -- Analyse
00314 // ---------------------------------------------------------------------------------//
00315 void TrackBuildingAnalyzer::analyze
00316 (
00317     const edm::Event& iEvent,
00318     const edm::EventSetup& iSetup,
00319     const TrackCandidate& candidate,
00320     const reco::BeamSpot& bs,
00321     const edm::ESHandle<MagneticField>& theMF,
00322     const edm::ESHandle<TransientTrackingRecHitBuilder>& theTTRHBuilder
00323 )
00324 {
00325 
00326     
00327     TSCBLBuilderNoMaterial tscblBuilder;
00328 
00329     //get parameters and errors from the candidate state
00330     TransientTrackingRecHit::RecHitPointer recHit = theTTRHBuilder->build(&*(candidate.recHits().second-1));
00331     TrajectoryStateOnSurface state = trajectoryStateTransform::transientState( candidate.trajectoryStateOnDet(), recHit->surface(), theMF.product());
00332     TrajectoryStateClosestToBeamLine tsAtClosestApproachTrackCand = tscblBuilder(*state.freeState(),bs);//as in TrackProducerAlgorithm
00333     if(!(tsAtClosestApproachTrackCand.isValid())) {
00334         edm::LogVerbatim("TrackBuilding") << "TrajectoryStateClosestToBeamLine not valid";
00335         return;
00336     }
00337     GlobalPoint  v0 = tsAtClosestApproachTrackCand.trackStateAtPCA().position();
00338     GlobalVector p = tsAtClosestApproachTrackCand.trackStateAtPCA().momentum();
00339     GlobalPoint  v(v0.x()-bs.x0(),v0.y()-bs.y0(),v0.z()-bs.z0());
00340 
00341     double pt           = sqrt(state.globalMomentum().perp2());
00342     double eta          = state.globalMomentum().eta();
00343     double phi          = state.globalMomentum().phi();
00344     double theta        = state.globalMomentum().theta();
00345     //double pm           = sqrt(state.globalMomentum().mag2());
00346     //double pz           = state.globalMomentum().z();
00347     //double qoverp       = tsAtClosestApproachTrackCand.trackStateAtPCA().charge()/p.mag();
00348     //double theta        = p.theta();
00349     //double lambda       = M_PI/2-p.theta();
00350     double numberOfHits = candidate.recHits().second-candidate.recHits().first;
00351     double dxy          = (-v.x()*sin(p.phi())+v.y()*cos(p.phi()));
00352     
00353     double dz           = v.z() - (v.x()*p.x()+v.y()*p.y())/p.perp() * p.z()/p.perp();
00354     
00355     if (doTCPlots){
00356       // fill the ME's
00357       TrackCandQ->Fill( state.charge() );
00358       TrackCandPt->Fill( pt );
00359       TrackCandEta->Fill( eta );
00360       TrackCandPhi->Fill( phi );
00361       TrackCandTheta->Fill( theta );
00362       TrackCandDxy->Fill( dxy );
00363       TrackCandDz->Fill( dz );
00364       NumberOfRecHitsPerTrackCand->Fill( numberOfHits );
00365       NumberOfRecHitsPerTrackCandVsEtaProfile->Fill( eta, numberOfHits );
00366       NumberOfRecHitsPerTrackCandVsPhiProfile->Fill( phi, numberOfHits );
00367     }
00368 }