Reference Manual

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    CosmicSplitterValidation
00004 // Class:      CosmicSplitterValidation
00005 // 
00013 //
00014 // Original Author:  Nhan Tran
00015 //         Created:  Mon Jul 16m 16:56:34 CDT 2007
00016 // $Id: CosmicSplitterValidation.cc,v 1.12 2011/12/20 15:11:41 mussgill Exp $
00017 //
00018 //
00019 
00020 // system include files
00021 #include "FWCore/Framework/interface/EDAnalyzer.h"
00022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00023 #include "FWCore/Framework/interface/MakerMacros.h"
00024 
00025 #include <algorithm>
00026 
00027 #include "FWCore/Framework/interface/ESHandle.h"
00028 #include "FWCore/Framework/interface/EventSetup.h"
00029 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00030 
00031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00032 #include "FWCore/Utilities/interface/InputTag.h"
00033 
00034 
00035 #include "FWCore/Framework/interface/EDProducer.h"
00036 #include "FWCore/Framework/interface/Event.h"
00037 #include "FWCore/Framework/interface/EventSetup.h"
00038 #include "FWCore/Framework/interface/ESHandle.h"
00039 #include "CommonTools/UtilAlgos/interface/TFileService.h"
00040 #include "FWCore/ServiceRegistry/interface/Service.h"
00041 
00042 
00043 #include "DataFormats/Common/interface/View.h"
00044 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00045 #include "DataFormats/TrackReco/interface/Track.h"
00046 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
00047 #include "DataFormats/TrackingRecHit/interface/InvalidTrackingRecHit.h"
00048 #include "DataFormats/MuonReco/interface/Muon.h"
00049 #include "DataFormats/MuonReco/interface/MuonFwd.h"
00050 
00051 #include "TrackingTools/TransientTrack/interface/BasicTransientTrack.h"
00052 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h" 
00053 #include "DataFormats/TrackReco/interface/Track.h"
00054 #include "DataFormats/TrackReco/interface/TrackFwd.h" 
00055 #include "DataFormats/Common/interface/RefToBase.h" 
00056 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHitBuilder.h"
00057 
00058 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h" 
00059 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00060 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
00061 #include "TrackingTools/TrajectoryState/interface/CopyUsingClone.h" 
00062 #include "RecoVertex/VertexTools/interface/PerigeeLinearizedTrackState.h" 
00063 
00064 #include <TTree.h>
00065 
00066 //
00067 // class decleration
00068 //
00069 
00070 class CosmicSplitterValidation : public edm::EDAnalyzer {
00071 public:
00072         explicit CosmicSplitterValidation(const edm::ParameterSet&);
00073         ~CosmicSplitterValidation();
00074         
00075         
00076 private:
00077         virtual void beginJob();
00078         virtual void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup);
00079         virtual void endJob() ;
00080         
00081         bool is_gold_muon(const edm::Event& e);
00082         
00083         edm::InputTag splitTracks_;
00084         edm::InputTag splitGlobalMuons_;
00085         edm::InputTag originalTracks_;
00086         edm::InputTag originalGlobalMuons_;
00087         
00088         bool checkIfGolden_;
00089         bool splitMuons_;
00090         int totalTracksToAnalyzer_;
00091         int goldenCtr;
00092         int twoTracksCtr;
00093         int goldenPlusTwoTracksCtr;
00094         int _passesTracksPlusMuonsCuts;
00095         
00096         edm::Service<TFileService> tfile;
00097         // ----------member data ---------------------------
00098         //std::vector<AlignTransform> m_align;
00099         // tree
00100         TTree* splitterTree_;
00101         // tree vars
00102         // split track variables
00103         double dcaX1_spl_, dcaY1_spl_, dcaZ1_spl_;
00104         double dcaX2_spl_, dcaY2_spl_, dcaZ2_spl_;
00105         double dxy1_spl_, dxy2_spl_, dz1_spl_, dz2_spl_;
00106         double theta1_spl_, theta2_spl_, phi1_spl_, phi2_spl_;
00107         double ddxy_spl_, ddz_spl_, dtheta_spl_, dphi_spl_;
00108         double pt1_spl_, pt2_spl_, dpt_spl_, p1_spl_, p2_spl_;
00109         double eta1_spl_, eta2_spl_, deta_spl_;
00110         double nHits1_spl_, nHitsPXB1_spl_, nHitsPXF1_spl_, nHitsTIB1_spl_;
00111         double nHitsTOB1_spl_, nHitsTID1_spl_, nHitsTEC1_spl_;
00112         double nHits2_spl_, nHitsPXB2_spl_, nHitsPXF2_spl_, nHitsTIB2_spl_;
00113         double nHitsTOB2_spl_, nHitsTID2_spl_, nHitsTEC2_spl_;
00114         // spl_it track errors
00115         double pt1Err_spl_, pt2Err_spl_;
00116         double theta1Err_spl_, theta2Err_spl_;
00117         double phi1Err_spl_, phi2Err_spl_;
00118         double d01Err_spl_, d02Err_spl_;
00119         double dz1Err_spl_, dz2Err_spl_;
00120         // original track variables
00121         double dcaX_org_, dcaY_org_, dcaZ_org_;
00122         double dxy_org_, dz_org_;
00123         double theta_org_, phi_org_, eta_org_, pt_org_, p_org_;
00124         double norchi2_org_;
00125         
00126         // split sta variables
00127         double dcaX1_sta_, dcaY1_sta_, dcaZ1_sta_;
00128         double dcaX2_sta_, dcaY2_sta_, dcaZ2_sta_;
00129         double dxy1_sta_, dxy2_sta_, dz1_sta_, dz2_sta_;
00130         double theta1_sta_, theta2_sta_, phi1_sta_, phi2_sta_;
00131         double ddxy_sta_, ddz_sta_, dtheta_sta_, dphi_sta_;
00132         double pt1_sta_, pt2_sta_, dpt_sta_, p1_sta_, p2_sta_;
00133         double eta1_sta_, eta2_sta_, deta_sta_;
00134         // split sta_ errors
00135         double pt1Err_sta_, pt2Err_sta_;
00136         double theta1Err_sta_, theta2Err_sta_;
00137         double phi1Err_sta_, phi2Err_sta_;
00138         double d01Err_sta_, d02Err_sta_;
00139         double dz1Err_sta_, dz2Err_sta_;
00140 
00141         // split glb_ variables
00142         double dcaX1_glb_, dcaY1_glb_, dcaZ1_glb_;
00143         double dcaX2_glb_, dcaY2_glb_, dcaZ2_glb_;
00144         double dxy1_glb_, dxy2_glb_, dz1_glb_, dz2_glb_;
00145         double theta1_glb_, theta2_glb_, phi1_glb_, phi2_glb_;
00146         double ddxy_glb_, ddz_glb_, dtheta_glb_, dphi_glb_;
00147         double pt1_glb_, pt2_glb_, dpt_glb_, p1_glb_, p2_glb_;
00148         double eta1_glb_, eta2_glb_, deta_glb_;
00149         double norchi1_glb_, norchi2_glb_;
00150         // split glb_ errors
00151         double pt1Err_glb_, pt2Err_glb_;
00152         double theta1Err_glb_, theta2Err_glb_;
00153         double phi1Err_glb_, phi2Err_glb_;
00154         double d01Err_glb_, d02Err_glb_;
00155         double dz1Err_glb_, dz2Err_glb_;
00156         // original glb muon variables
00157         double dcaX_orm_, dcaY_orm_, dcaZ_orm_;
00158         double dxy_orm_, dz_orm_;
00159         double theta_orm_, phi_orm_, eta_orm_, pt_orm_, p_orm_;
00160         double norchi2_orm_;
00161         
00162 };
00163 
00164 //
00165 // constants, enums and typedefs
00166 //
00167 
00168 //
00169 // static data member definitions
00170 //
00171 
00172 //
00173 // constructors and destructor
00174 //
00175 CosmicSplitterValidation::CosmicSplitterValidation(const edm::ParameterSet& iConfig):
00176         splitTracks_(iConfig.getParameter<edm::InputTag>("splitTracks")),
00177         splitGlobalMuons_(iConfig.getParameter<edm::InputTag>("splitGlobalMuons")),
00178         originalTracks_(iConfig.getParameter<edm::InputTag>("originalTracks")),
00179         originalGlobalMuons_(iConfig.getParameter<edm::InputTag>("originalGlobalMuons")),
00180         checkIfGolden_(iConfig.getParameter<bool>("checkIfGolden")),
00181         splitMuons_(iConfig.getParameter<bool> ("ifSplitMuons")),
00182         totalTracksToAnalyzer_(0),
00183         goldenCtr(0),
00184         twoTracksCtr(0),
00185         goldenPlusTwoTracksCtr(0),
00186         _passesTracksPlusMuonsCuts(0),
00187         splitterTree_(0),
00188         dcaX1_spl_(0), dcaY1_spl_(0), dcaZ1_spl_(0),
00189         dcaX2_spl_(0), dcaY2_spl_(0), dcaZ2_spl_(0),
00190         dxy1_spl_(0), dxy2_spl_(0), dz1_spl_(0), dz2_spl_(0),
00191         theta1_spl_(0), theta2_spl_(0), phi1_spl_(0), phi2_spl_(0),
00192         ddxy_spl_(0), ddz_spl_(0), dtheta_spl_(0), dphi_spl_(0),
00193         pt1_spl_(0), pt2_spl_(0), dpt_spl_(0), p1_spl_(0), p2_spl_(0),
00194         eta1_spl_(0), eta2_spl_(0), deta_spl_(0),
00195         nHits1_spl_(0), nHitsPXB1_spl_(0), nHitsPXF1_spl_(0), nHitsTIB1_spl_(0),
00196         nHitsTOB1_spl_(0), nHitsTID1_spl_(0), nHitsTEC1_spl_(0),
00197         nHits2_spl_(0), nHitsPXB2_spl_(0), nHitsPXF2_spl_(0), nHitsTIB2_spl_(0),
00198         nHitsTOB2_spl_(0), nHitsTID2_spl_(0), nHitsTEC2_spl_(0),
00199         pt1Err_spl_(0), pt2Err_spl_(0),
00200         theta1Err_spl_(0), theta2Err_spl_(0),
00201         phi1Err_spl_(0), phi2Err_spl_(0),
00202         d01Err_spl_(0), d02Err_spl_(0),
00203         dz1Err_spl_(0), dz2Err_spl_(0),
00204         dcaX_org_(0), dcaY_org_(0), dcaZ_org_(0),
00205         dxy_org_(0), dz_org_(0),
00206         theta_org_(0), phi_org_(0), eta_org_(0), pt_org_(0), p_org_(0),
00207         norchi2_org_(0),
00208         dcaX1_sta_(0), dcaY1_sta_(0), dcaZ1_sta_(0),
00209         dcaX2_sta_(0), dcaY2_sta_(0), dcaZ2_sta_(0),
00210         dxy1_sta_(0), dxy2_sta_(0), dz1_sta_(0), dz2_sta_(0),
00211         theta1_sta_(0), theta2_sta_(0), phi1_sta_(0), phi2_sta_(0),
00212         ddxy_sta_(0), ddz_sta_(0), dtheta_sta_(0), dphi_sta_(0),
00213         pt1_sta_(0), pt2_sta_(0), dpt_sta_(0), p1_sta_(0), p2_sta_(0),
00214         eta1_sta_(0), eta2_sta_(0), deta_sta_(0),
00215         pt1Err_sta_(0), pt2Err_sta_(0),
00216         theta1Err_sta_(0), theta2Err_sta_(0),
00217         phi1Err_sta_(0), phi2Err_sta_(0),
00218         d01Err_sta_(0), d02Err_sta_(0),
00219         dz1Err_sta_(0), dz2Err_sta_(0),
00220         dcaX1_glb_(0), dcaY1_glb_(0), dcaZ1_glb_(0),
00221         dcaX2_glb_(0), dcaY2_glb_(0), dcaZ2_glb_(0),
00222         dxy1_glb_(0), dxy2_glb_(0), dz1_glb_(0), dz2_glb_(0),
00223         theta1_glb_(0), theta2_glb_(0), phi1_glb_(0), phi2_glb_(0),
00224         ddxy_glb_(0), ddz_glb_(0), dtheta_glb_(0), dphi_glb_(0),
00225         pt1_glb_(0), pt2_glb_(0), dpt_glb_(0), p1_glb_(0), p2_glb_(0),
00226         eta1_glb_(0), eta2_glb_(0), deta_glb_(0),
00227         norchi1_glb_(0), norchi2_glb_(0),
00228         pt1Err_glb_(0), pt2Err_glb_(0),
00229         theta1Err_glb_(0), theta2Err_glb_(0),
00230         phi1Err_glb_(0), phi2Err_glb_(0),
00231         d01Err_glb_(0), d02Err_glb_(0),
00232         dz1Err_glb_(0), dz2Err_glb_(0),
00233         dcaX_orm_(0), dcaY_orm_(0), dcaZ_orm_(0),
00234         dxy_orm_(0), dz_orm_(0),
00235         theta_orm_(0), phi_orm_(0), eta_orm_(0), pt_orm_(0), p_orm_(0),
00236         norchi2_orm_(0)
00237 {
00238         
00239 }
00240 
00241 
00242 CosmicSplitterValidation::~CosmicSplitterValidation()
00243 {}
00244 
00245 
00246 //
00247 // member functions
00248 //
00249 
00250 // ------------ method called to for each event  ------------
00251 void CosmicSplitterValidation::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00252 {
00253         
00254         // check if golden muon
00255         bool isGolden = true;
00256         if (checkIfGolden_) isGolden = is_gold_muon( iEvent );
00257 
00258         // grab collections
00259         edm::Handle<std::vector<reco::Track> > tracks;
00260         edm::Handle<reco::MuonCollection> globalMuons;
00261         edm::Handle<reco::MuonCollection> originalGlobalMuons;
00262         edm::Handle<std::vector<reco::Track> > originalTracks;
00263         iEvent.getByLabel(splitTracks_, tracks);
00264         iEvent.getByLabel(originalTracks_, originalTracks);
00265         if (splitMuons_){
00266                 iEvent.getByLabel(splitGlobalMuons_, globalMuons);
00267                 iEvent.getByLabel(originalGlobalMuons_, originalGlobalMuons);
00268         }
00269 
00270         const int kBPIX = PixelSubdetector::PixelBarrel;
00271         const int kFPIX = PixelSubdetector::PixelEndcap;
00272 
00273         
00274         totalTracksToAnalyzer_ = totalTracksToAnalyzer_ + tracks->size();
00275         if (isGolden) goldenCtr++;
00276         if (tracks->size() == 2) twoTracksCtr++;
00277         if (tracks->size() == 2 && originalTracks->size() == 1 && isGolden){
00278                 goldenPlusTwoTracksCtr++;
00279                 
00280                 int gmCtr = 0; 
00281                 bool topGlobalMuonFlag = false;
00282                 bool bottomGlobalMuonFlag = false;
00283                 int topGlobalMuon = -1;
00284                 int bottomGlobalMuon = -1;
00285                 double topGlobalMuonNorchi2 = 1e10;
00286                 double bottomGlobalMuonNorchi2 = 1e10;
00287                 
00288                 if (splitMuons_){
00289                         // check if split global muons are good
00290                         for (std::vector<reco::Muon>::const_iterator gmI = globalMuons->begin(); gmI != globalMuons->end(); gmI++){
00291                                 
00292                                 if ( gmI->isTrackerMuon() && gmI->isStandAloneMuon() && gmI->isGlobalMuon() ){
00293                                         
00294                                         reco::TrackRef trackerTrackRef1( tracks, 0 );
00295                                         reco::TrackRef trackerTrackRef2( tracks, 1 );
00296                                         
00297                                         if (gmI->innerTrack() == trackerTrackRef1){
00298                                                 if (gmI->globalTrack()->normalizedChi2() < topGlobalMuonNorchi2){
00299                                                         topGlobalMuonFlag = true;
00300                                                         topGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
00301                                                         topGlobalMuon = gmCtr;
00302                                                 }
00303                                         }
00304                                         if (gmI->innerTrack() == trackerTrackRef2){
00305                                                 if (gmI->globalTrack()->normalizedChi2() < bottomGlobalMuonNorchi2){
00306                                                         bottomGlobalMuonFlag = true;
00307                                                         bottomGlobalMuonNorchi2 = gmI->globalTrack()->normalizedChi2();
00308                                                         bottomGlobalMuon = gmCtr;
00309                                                 }
00310                                         }
00311                                 }
00312                                 gmCtr++;
00313                         }                       
00314                 }
00315                 
00316                 if ( (!splitMuons_) || (splitMuons_ && bottomGlobalMuonFlag && topGlobalMuonFlag) ){
00317                         
00318                         _passesTracksPlusMuonsCuts++;
00319                         
00320                         // split tracks calculations
00321                         reco::Track track1 = tracks->at(0);
00322                         reco::Track track2 = tracks->at(1);
00323                         
00324                         math::XYZPoint dca1 = track1.referencePoint();
00325                         math::XYZPoint dca2 = track2.referencePoint();
00326                         
00327                         // looping through the hits for track 1
00328                         double Nrechits1 =0;
00329                         double nhitinTIB1 =0; 
00330                         double nhitinTOB1 =0; 
00331                         double nhitinTID1 =0; 
00332                         double nhitinTEC1 =0; 
00333                         double nhitinBPIX1 =0;
00334                         double nhitinFPIX1 =0;    
00335                         for (trackingRecHit_iterator iHit = track1.recHitsBegin(); iHit != track1.recHitsEnd(); ++iHit) {
00336                                 
00337                                 if((*iHit)->isValid()) {       
00338                                         
00339                                         Nrechits1++;
00340                                         
00341                                         int type =(*iHit)->geographicalId().subdetId();
00342                                         
00343                                         if(type==int(StripSubdetector::TIB)){++nhitinTIB1;}
00344                                         if(type==int(StripSubdetector::TOB)){++nhitinTOB1;}
00345                                         if(type==int(StripSubdetector::TID)){++nhitinTID1;}
00346                                         if(type==int(StripSubdetector::TEC)){++nhitinTEC1;}
00347                                         if(type==int(                kBPIX)){++nhitinBPIX1;}
00348                                         if(type==int(                kFPIX)){++nhitinFPIX1;}
00349                                         
00350                                 }
00351                         }    
00352                         nHits1_spl_ = Nrechits1;
00353                         nHitsTIB1_spl_ = nhitinTIB1; 
00354                         nHitsTOB1_spl_ = nhitinTOB1; 
00355                         nHitsTID1_spl_ = nhitinTID1; 
00356                         nHitsTEC1_spl_ = nhitinTEC1; 
00357                         nHitsPXB1_spl_ = nhitinBPIX1;
00358                         nHitsPXF1_spl_ = nhitinFPIX1;  
00359                         
00360                         // looping through the hits for track 2
00361                         double Nrechits2 =0;
00362                         double nhitinTIB2 =0; 
00363                         double nhitinTOB2 =0; 
00364                         double nhitinTID2 =0; 
00365                         double nhitinTEC2 =0; 
00366                         double nhitinBPIX2 =0;
00367                         double nhitinFPIX2 =0;    
00368                         for (trackingRecHit_iterator iHit = track2.recHitsBegin(); iHit != track2.recHitsEnd(); ++iHit) {
00369                                 
00370                                 if((*iHit)->isValid()) {       
00371                                         
00372                                         Nrechits2++;
00373                                         
00374                                         int type =(*iHit)->geographicalId().subdetId();
00375                                         
00376                                         if(type==int(StripSubdetector::TIB)){++nhitinTIB2;}
00377                                         if(type==int(StripSubdetector::TOB)){++nhitinTOB2;}
00378                                         if(type==int(StripSubdetector::TID)){++nhitinTID2;}
00379                                         //\\if(type==int(StripSubdetector::TEC)){++nhitinTEC2;}
00380                                         if(type==int(                kBPIX)){++nhitinBPIX2;}
00381                                         if(type==int(                kFPIX)){++nhitinFPIX2;}
00382                                         
00383                                 }
00384                         }    
00385                         nHits2_spl_ = Nrechits2;
00386                         nHitsTIB2_spl_ = nhitinTIB2; 
00387                         nHitsTOB2_spl_ = nhitinTOB2; 
00388                         nHitsTID2_spl_ = nhitinTID2; 
00389                         nHitsTEC2_spl_ = nhitinTEC2; 
00390                         nHitsPXB2_spl_ = nhitinBPIX2;
00391                         nHitsPXF2_spl_ = nhitinFPIX2;    
00392                         
00393                         
00394                         double dtheta_Val = track1.theta() - track2.theta();
00395                         double dphi_Val = track1.phi() - track2.phi();
00396                         double ddxy_Val = track1.d0() - track2.d0();
00397                         double ddz_Val = track1.dz() - track2.dz();
00398                         double dpt_Val = track1.pt() - track2.pt();
00399                         
00400                         // original tracks calculations
00401                         reco::Track origTrack = originalTracks->at(0);
00402                         math::XYZPoint dca_org = origTrack.referencePoint();
00403                         
00404                         // fill tree
00405                         // split tracks
00406                         dcaX1_spl_ = dca1.x(); 
00407                         dcaY1_spl_ = dca1.y();
00408                         dcaZ1_spl_ = dca1.z();
00409                         dcaX2_spl_ = dca2.x(); 
00410                         dcaY2_spl_ = dca2.y();
00411                         dcaZ2_spl_ = dca2.z();
00412                         dxy1_spl_ = track1.d0();
00413                         dxy2_spl_ = track2.d0();
00414                         dz1_spl_ = track1.dz();
00415                         dz2_spl_ = track2.dz();
00416                         d01Err_spl_ = track1.d0Error();
00417                         d02Err_spl_ = track2.d0Error();
00418                         dz1Err_spl_ = track1.dzError();
00419                         dz2Err_spl_ = track2.dzError();
00420                         theta1_spl_ = track1.theta();
00421                         theta2_spl_ = track2.theta();
00422                         theta1Err_spl_ = track1.thetaError();
00423                         theta2Err_spl_ = track2.thetaError();
00424                         phi1_spl_ = track1.phi();
00425                         phi2_spl_ = track2.phi();
00426                         phi1Err_spl_ = track1.phiError();
00427                         phi2Err_spl_ = track2.phiError();
00428                         ddxy_spl_ = ddxy_Val;
00429                         ddz_spl_ = ddz_Val;
00430                         dtheta_spl_ = dtheta_Val;
00431                         dphi_spl_ = dphi_Val;
00432                         pt1_spl_ = track1.pt();
00433                         pt2_spl_ = track2.pt();
00434                         p1_spl_ = track1.p();
00435                         p2_spl_ = track2.p();
00436                         pt1Err_spl_ = track1.ptError();
00437                         pt2Err_spl_ = track2.ptError();
00438                         dpt_spl_ = dpt_Val;
00439                         eta1_spl_ = track1.eta();
00440                         eta2_spl_ = track2.eta();
00441                         deta_spl_ = eta1_spl_ - eta2_spl_;
00442                         
00443                         // original tracks
00444                         dcaX_org_ = dca_org.x();
00445                         dcaY_org_ = dca_org.y();
00446                         dcaZ_org_ = dca_org.z();
00447                         dxy_org_ = origTrack.d0();
00448                         dz_org_ = origTrack.dz();
00449                         theta_org_ = origTrack.theta();
00450                         phi_org_ = origTrack.phi();
00451                         eta_org_ = origTrack.eta();
00452                         pt_org_ = origTrack.pt();
00453                         p_org_ = origTrack.p();
00454                         norchi2_org_ = origTrack.normalizedChi2();
00455                         
00456                         // split muons calculations
00457                         if (splitMuons_){
00458                                 
00459                                 reco::Muon muonTop = globalMuons->at( topGlobalMuon );
00460                                 reco::Muon muonBottom = globalMuons->at( bottomGlobalMuon );                            
00461                                 
00462                                 reco::TrackRef glb1 = muonTop.globalTrack();
00463                                 reco::TrackRef glb2 = muonBottom.globalTrack();
00464                                 reco::TrackRef sta1 = muonTop.outerTrack();
00465                                 reco::TrackRef sta2 = muonBottom.outerTrack();
00466                                 
00467                                 // standalone muon variables
00468                                 dcaX1_sta_ = sta1->referencePoint().x();
00469                                 dcaY1_sta_ = sta1->referencePoint().y();
00470                                 dcaZ1_sta_ = sta1->referencePoint().z();
00471                                 dcaX2_sta_ = sta2->referencePoint().x();
00472                                 dcaY2_sta_ = sta2->referencePoint().y();
00473                                 dcaZ2_sta_ = sta2->referencePoint().z();
00474                                 dxy1_sta_ = sta1->d0();
00475                                 dxy2_sta_ = sta2->d0();
00476                                 dz1_sta_ = sta1->dz();
00477                                 dz2_sta_ = sta2->dz();
00478                                 d01Err_sta_ = sta1->d0Error();
00479                                 d02Err_sta_ = sta2->d0Error();
00480                                 dz1Err_sta_ = sta1->dzError();
00481                                 dz2Err_sta_ = sta2->dzError();
00482                                 theta1_sta_ = sta1->theta();
00483                                 theta2_sta_ = sta2->theta();
00484                                 theta1Err_sta_ = sta1->thetaError();
00485                                 theta2Err_sta_ = sta2->thetaError();
00486                                 phi1_sta_ = sta1->phi();
00487                                 phi2_sta_ = sta2->phi();
00488                                 phi1Err_sta_ = sta1->phiError();
00489                                 phi2Err_sta_ = sta2->phiError();
00490                                 ddxy_sta_ = sta1->d0() - sta2->d0();
00491                                 ddz_sta_ = sta1->dz() - sta2->dz();
00492                                 dtheta_sta_ = sta1->theta() - sta2->theta();
00493                                 dphi_sta_ = sta1->phi() - sta2->phi();
00494                                 pt1_sta_ = sta1->pt();
00495                                 pt2_sta_ = sta2->pt();
00496                                 p1_sta_ = sta1->p();
00497                                 p2_sta_ = sta2->p();
00498                                 pt1Err_sta_ = sta1->ptError();
00499                                 pt2Err_sta_ = sta2->ptError();
00500                                 dpt_sta_ = sta1->pt() - sta2->pt();
00501                                 eta1_sta_ = sta1->eta();
00502                                 eta2_sta_ = sta2->eta();
00503                                 deta_sta_ = eta1_sta_ - eta2_sta_;
00504                                 
00505                                 // global muon variables
00506                                 dcaX1_glb_ = glb1->referencePoint().x();
00507                                 dcaY1_glb_ = glb1->referencePoint().y();
00508                                 dcaZ1_glb_ = glb1->referencePoint().z();
00509                                 dcaX2_glb_ = glb2->referencePoint().x();
00510                                 dcaY2_glb_ = glb2->referencePoint().y();
00511                                 dcaZ2_glb_ = glb2->referencePoint().z();
00512                                 dxy1_glb_ = glb1->d0();
00513                                 dxy2_glb_ = glb2->d0();
00514                                 dz1_glb_ = glb1->dz();
00515                                 dz2_glb_ = glb2->dz();
00516                                 d01Err_glb_ = glb1->d0Error();
00517                                 d02Err_glb_ = glb2->d0Error();
00518                                 dz1Err_glb_ = glb1->dzError();
00519                                 dz2Err_glb_ = glb2->dzError();
00520                                 theta1_glb_ = glb1->theta();
00521                                 theta2_glb_ = glb2->theta();
00522                                 theta1Err_glb_ = glb1->thetaError();
00523                                 theta2Err_glb_ = glb2->thetaError();
00524                                 phi1_glb_ = glb1->phi();
00525                                 phi2_glb_ = glb2->phi();
00526                                 phi1Err_glb_ = glb1->phiError();
00527                                 phi2Err_glb_ = glb2->phiError();
00528                                 ddxy_glb_ = glb1->d0() - glb2->d0();
00529                                 ddz_glb_ = glb1->dz() - glb2->dz();
00530                                 dtheta_glb_ = glb1->theta() - glb2->theta();
00531                                 dphi_glb_ = glb1->phi() - glb2->phi();
00532                                 pt1_glb_ = glb1->pt();
00533                                 pt2_glb_ = glb2->pt();
00534                                 p1_glb_ = glb1->p();
00535                                 p2_glb_ = glb2->p();
00536                                 pt1Err_glb_ = glb1->ptError();
00537                                 pt2Err_glb_ = glb2->ptError();
00538                                 dpt_glb_ = glb1->pt() - glb2->pt();
00539                                 eta1_glb_ = glb1->eta();
00540                                 eta2_glb_ = glb2->eta();
00541                                 deta_glb_ = eta1_glb_ - eta2_glb_;
00542                                 norchi1_glb_ = glb1->normalizedChi2();
00543                                 norchi2_glb_ = glb2->normalizedChi2();
00544                                 
00545                         }
00546                         
00547                         
00548                         splitterTree_->Fill();
00549                 }
00550         }
00551         
00552         
00553 }
00554 
00555 
00556 // ------------ method called once each job just before starting event loop  ------------
00557 void CosmicSplitterValidation::beginJob()
00558 {
00559         edm::LogInfo("beginJob") << "Begin Job" << std::endl;
00560         
00561         splitterTree_ = tfile->make<TTree>("splitterTree","splitterTree");
00562         
00563         // split track variables
00564         splitterTree_->Branch("dcaX1_spl", &dcaX1_spl_, "dcaX1_spl/D");
00565         splitterTree_->Branch("dcaY1_spl", &dcaY1_spl_, "dcaY1_spl/D");
00566         splitterTree_->Branch("dcaZ1_spl", &dcaZ1_spl_, "dcaZ1_spl/D");
00567         splitterTree_->Branch("dcaX2_spl", &dcaX2_spl_, "dcaX2_spl/D");
00568         splitterTree_->Branch("dcaY2_spl", &dcaY2_spl_, "dcaY2_spl/D");
00569         splitterTree_->Branch("dcaZ2_spl", &dcaZ2_spl_, "dcaZ2_spl/D");
00570         splitterTree_->Branch("dxy1_spl", &dxy1_spl_, "dxy1_spl/D");
00571         splitterTree_->Branch("dz1_spl", &dz1_spl_, "dz1_spl/D");
00572         splitterTree_->Branch("dxy2_spl", &dxy2_spl_, "dxy2_spl/D");
00573         splitterTree_->Branch("dz2_spl", &dz2_spl_, "dz2_spl/D");
00574         splitterTree_->Branch("theta1_spl", &theta1_spl_, "theta1_spl/D");
00575         splitterTree_->Branch("theta2_spl", &theta2_spl_, "theta2_spl/D");
00576         splitterTree_->Branch("phi1_spl", &phi1_spl_, "phi1_spl/D");
00577         splitterTree_->Branch("phi2_spl", &phi2_spl_, "phi2_spl/D");
00578         splitterTree_->Branch("ddxy_spl", &ddxy_spl_, "ddxy_spl/D");
00579         splitterTree_->Branch("ddz_spl", &ddz_spl_, "ddz_spl/D");
00580         splitterTree_->Branch("dphi_spl", &dphi_spl_, "dphi_spl/D");
00581         splitterTree_->Branch("dtheta_spl", &dtheta_spl_, "dtheta_spl/D");
00582         splitterTree_->Branch("pt1_spl", &pt1_spl_, "pt1_spl/D");
00583         splitterTree_->Branch("pt2_spl", &pt2_spl_, "pt2_spl/D");
00584         splitterTree_->Branch("dpt_spl", &dpt_spl_, "dpt_spl/D");
00585         splitterTree_->Branch("p1_spl", &p1_spl_, "p1_spl/D");
00586         splitterTree_->Branch("p2_spl", &p2_spl_, "p2_spl/D");
00587         splitterTree_->Branch("eta1_spl", &eta1_spl_, "eta1_spl/D");
00588         splitterTree_->Branch("eta2_spl", &eta2_spl_, "eta2_spl/D");
00589         splitterTree_->Branch("deta_spl", &deta_spl_, "deta_spl/D");
00590         splitterTree_->Branch("nHits1_spl", &nHits1_spl_, "nHits1_spl/D");
00591         splitterTree_->Branch("nHitsPXB1_spl", &nHitsPXB1_spl_, "nHitsPXB1_spl/D");
00592         splitterTree_->Branch("nHitsPXF1_spl", &nHitsPXF1_spl_, "nHitsPXF1_spl/D");
00593         splitterTree_->Branch("nHitsTIB1_spl", &nHitsTIB1_spl_, "nHitsTIB1_spl/D");
00594         splitterTree_->Branch("nHitsTOB1_spl", &nHitsTOB1_spl_, "nHitsTOB1_spl/D");
00595         splitterTree_->Branch("nHitsTID1_spl", &nHitsTID1_spl_, "nHitsTID1_spl/D");
00596         splitterTree_->Branch("nHitsTEC1_spl", &nHitsTEC1_spl_, "nHitsTEC1_spl/D");
00597         splitterTree_->Branch("nHits2_spl", &nHits2_spl_, "nHits2_spl/D");
00598         splitterTree_->Branch("nHitsPXB2_spl", &nHitsPXB2_spl_, "nHitsPXB2_spl/D");
00599         splitterTree_->Branch("nHitsPXF2_spl", &nHitsPXF2_spl_, "nHitsPXF2_spl/D");
00600         splitterTree_->Branch("nHitsTIB2_spl", &nHitsTIB2_spl_, "nHitsTIB2_spl/D");
00601         splitterTree_->Branch("nHitsTOB2_spl", &nHitsTOB2_spl_, "nHitsTOB2_spl/D");
00602         splitterTree_->Branch("nHitsTID2_spl", &nHitsTID2_spl_, "nHitsTID2_spl/D");
00603         splitterTree_->Branch("nHitsTEC2_spl", &nHitsTEC2_spl_, "nHitsTEC2_spl/D");
00604         
00605         
00606         splitterTree_->Branch("d01Err_spl", &d01Err_spl_, "d01Err_spl/D");
00607         splitterTree_->Branch("d02Err_spl", &d02Err_spl_, "d02Err_spl/D");
00608         splitterTree_->Branch("dz1Err_spl", &dz1Err_spl_, "dz1Err_spl/D");
00609         splitterTree_->Branch("dz2Err_spl", &dz2Err_spl_, "dz2Err_spl/D");
00610         splitterTree_->Branch("phi1Err_spl", &phi1Err_spl_, "phi1Err_spl/D");
00611         splitterTree_->Branch("phi2Err_spl", &phi2Err_spl_, "phi2Err_spl/D");
00612         splitterTree_->Branch("theta1Err_spl", &theta1Err_spl_, "theta1Err_spl/D");
00613         splitterTree_->Branch("theta2Err_spl", &theta2Err_spl_, "theta2Err_spl/D");
00614         splitterTree_->Branch("pt1Err_spl", &pt1Err_spl_, "pt1Err_spl/D");
00615         splitterTree_->Branch("pt2Err_spl", &pt2Err_spl_, "pt2Err_spl/D");
00616         
00617         splitterTree_->Branch("dcaX_org", &dcaX_org_, "dcaX_org/D");
00618         splitterTree_->Branch("dcaY_org", &dcaY_org_, "dcaY_org/D");
00619         splitterTree_->Branch("dcaZ_org", &dcaZ_org_, "dcaZ_org/D");
00620         splitterTree_->Branch("dxy_org", &dxy_org_, "dxy_org/D");
00621         splitterTree_->Branch("dz_org", &dz_org_, "dz_org/D");
00622         splitterTree_->Branch("theta_org", &theta_org_, "theta_org/D");
00623         splitterTree_->Branch("phi_org", &phi_org_, "phi_org/D");
00624         splitterTree_->Branch("eta_org", &eta_org_, "eta_org/D");
00625         splitterTree_->Branch("pt_org", &pt_org_, "pt_org/D");
00626         splitterTree_->Branch("p_org", &p_org_, "p_org/D");
00627         splitterTree_->Branch("norchi2_org", &norchi2_org_, "norchi2_org/D");
00628         
00629         if (splitMuons_){
00630                 
00631                 // standalone split 
00632                 splitterTree_->Branch("dcaX1_sta", &dcaX1_sta_, "dcaX1_sta/D");
00633                 splitterTree_->Branch("dcaY1_sta", &dcaY1_sta_, "dcaY1_sta/D");
00634                 splitterTree_->Branch("dcaZ1_sta", &dcaZ1_sta_, "dcaZ1_sta/D");
00635                 splitterTree_->Branch("dcaX2_sta", &dcaX2_sta_, "dcaX2_sta/D");
00636                 splitterTree_->Branch("dcaY2_sta", &dcaY2_sta_, "dcaY2_sta/D");
00637                 splitterTree_->Branch("dcaZ2_sta", &dcaZ2_sta_, "dcaZ2_sta/D");
00638                 splitterTree_->Branch("dxy1_sta", &dxy1_sta_, "dxy1_sta/D");
00639                 splitterTree_->Branch("dz1_sta", &dz1_sta_, "dz1_sta/D");
00640                 splitterTree_->Branch("dxy2_sta", &dxy2_sta_, "dxy2_sta/D");
00641                 splitterTree_->Branch("dz2_sta", &dz2_sta_, "dz2_sta/D");               
00642                 splitterTree_->Branch("theta1_sta", &theta1_sta_, "theta1_sta/D");
00643                 splitterTree_->Branch("theta2_sta", &theta2_sta_, "theta2_sta/D");
00644                 splitterTree_->Branch("phi1_sta", &phi1_sta_, "phi1_sta/D");
00645                 splitterTree_->Branch("phi2_sta", &phi2_sta_, "phi2_sta/D");
00646                 splitterTree_->Branch("ddxy_sta", &ddxy_sta_, "ddxy_sta/D");
00647                 splitterTree_->Branch("ddz_sta", &ddz_sta_, "ddz_sta/D");
00648                 splitterTree_->Branch("dphi_sta", &dphi_sta_, "dphi_sta/D");
00649                 splitterTree_->Branch("dtheta_sta", &dtheta_sta_, "dtheta_sta/D");
00650                 splitterTree_->Branch("pt1_sta", &pt1_sta_, "pt1_sta/D");
00651                 splitterTree_->Branch("pt2_sta", &pt2_sta_, "pt2_sta/D");
00652                 splitterTree_->Branch("dpt_sta", &dpt_sta_, "dpt_sta/D");
00653                 splitterTree_->Branch("p1_sta", &p1_sta_, "p1_sta/D");
00654                 splitterTree_->Branch("p2_sta", &p2_sta_, "p2_sta/D");
00655                 splitterTree_->Branch("eta1_sta", &eta1_sta_, "eta1_sta/D");
00656                 splitterTree_->Branch("eta2_sta", &eta2_sta_, "eta2_sta/D");
00657                 splitterTree_->Branch("deta_sta", &deta_sta_, "deta_sta/D");
00658                 
00659                 splitterTree_->Branch("d01Err_sta", &d01Err_sta_, "d01Err_sta/D");
00660                 splitterTree_->Branch("d02Err_sta", &d02Err_sta_, "d02Err_sta/D");
00661                 splitterTree_->Branch("dz1Err_sta", &dz1Err_sta_, "dz1Err_sta/D");
00662                 splitterTree_->Branch("dz2Err_sta", &dz2Err_sta_, "dz2Err_sta/D");
00663                 splitterTree_->Branch("phi1Err_sta", &phi1Err_sta_, "phi1Err_sta/D");
00664                 splitterTree_->Branch("phi2Err_sta", &phi2Err_sta_, "phi2Err_sta/D");
00665                 splitterTree_->Branch("theta1Err_sta", &theta1Err_sta_, "theta1Err_sta/D");
00666                 splitterTree_->Branch("theta2Err_sta", &theta2Err_sta_, "theta2Err_sta/D");
00667                 splitterTree_->Branch("pt1Err_sta", &pt1Err_sta_, "pt1Err_sta/D");
00668                 splitterTree_->Branch("pt2Err_sta", &pt2Err_sta_, "pt2Err_sta/D");
00669                 
00670                 // global split 
00671                 splitterTree_->Branch("dcaX1_glb", &dcaX1_glb_, "dcaX1_glb/D");
00672                 splitterTree_->Branch("dcaY1_glb", &dcaY1_glb_, "dcaY1_glb/D");
00673                 splitterTree_->Branch("dcaZ1_glb", &dcaZ1_glb_, "dcaZ1_glb/D");
00674                 splitterTree_->Branch("dcaX2_glb", &dcaX2_glb_, "dcaX2_glb/D");
00675                 splitterTree_->Branch("dcaY2_glb", &dcaY2_glb_, "dcaY2_glb/D");
00676                 splitterTree_->Branch("dcaZ2_glb", &dcaZ2_glb_, "dcaZ2_glb/D");
00677                 splitterTree_->Branch("dxy1_glb", &dxy1_glb_, "dxy1_glb/D");
00678                 splitterTree_->Branch("dz1_glb", &dz1_glb_, "dz1_glb/D");
00679                 splitterTree_->Branch("dxy2_glb", &dxy2_glb_, "dxy2_glb/D");
00680                 splitterTree_->Branch("dz2_glb", &dz2_glb_, "dz2_glb/D");               
00681                 splitterTree_->Branch("theta1_glb", &theta1_glb_, "theta1_glb/D");
00682                 splitterTree_->Branch("theta2_glb", &theta2_glb_, "theta2_glb/D");
00683                 splitterTree_->Branch("phi1_glb", &phi1_glb_, "phi1_glb/D");
00684                 splitterTree_->Branch("phi2_glb", &phi2_glb_, "phi2_glb/D");
00685                 splitterTree_->Branch("ddxy_glb", &ddxy_glb_, "ddxy_glb/D");
00686                 splitterTree_->Branch("ddz_glb", &ddz_glb_, "ddz_glb/D");
00687                 splitterTree_->Branch("dphi_glb", &dphi_glb_, "dphi_glb/D");
00688                 splitterTree_->Branch("dtheta_glb", &dtheta_glb_, "dtheta_glb/D");
00689                 splitterTree_->Branch("pt1_glb", &pt1_glb_, "pt1_glb/D");
00690                 splitterTree_->Branch("pt2_glb", &pt2_glb_, "pt2_glb/D");
00691                 splitterTree_->Branch("dpt_glb", &dpt_glb_, "dpt_glb/D");
00692                 splitterTree_->Branch("p1_glb", &p1_glb_, "p1_glb/D");
00693                 splitterTree_->Branch("p2_glb", &p2_glb_, "p2_glb/D");
00694                 splitterTree_->Branch("eta1_glb", &eta1_glb_, "eta1_glb/D");
00695                 splitterTree_->Branch("eta2_glb", &eta2_glb_, "eta2_glb/D");
00696                 splitterTree_->Branch("deta_glb", &deta_glb_, "deta_glb/D");
00697                 splitterTree_->Branch("norchi1_glb", &norchi1_glb_, "norchi1_glb/D");
00698                 splitterTree_->Branch("norchi2_glb", &norchi2_glb_, "norchi2_glb/D");
00699 
00700                 splitterTree_->Branch("d01Err_glb", &d01Err_glb_, "d01Err_glb/D");
00701                 splitterTree_->Branch("d02Err_glb", &d02Err_glb_, "d02Err_glb/D");
00702                 splitterTree_->Branch("dz1Err_glb", &dz1Err_glb_, "dz1Err_glb/D");
00703                 splitterTree_->Branch("dz2Err_glb", &dz2Err_glb_, "dz2Err_glb/D");
00704                 splitterTree_->Branch("phi1Err_glb", &phi1Err_glb_, "phi1Err_glb/D");
00705                 splitterTree_->Branch("phi2Err_glb", &phi2Err_glb_, "phi2Err_glb/D");
00706                 splitterTree_->Branch("theta1Err_glb", &theta1Err_glb_, "theta1Err_glb/D");
00707                 splitterTree_->Branch("theta2Err_glb", &theta2Err_glb_, "theta2Err_glb/D");
00708                 splitterTree_->Branch("pt1Err_glb", &pt1Err_glb_, "pt1Err_glb/D");
00709                 splitterTree_->Branch("pt2Err_glb", &pt2Err_glb_, "pt2Err_glb/D");
00710 
00711         }
00712         
00713         totalTracksToAnalyzer_ = 0;
00714         goldenCtr = 0;
00715         twoTracksCtr = 0;
00716         goldenPlusTwoTracksCtr = 0;
00717         _passesTracksPlusMuonsCuts = 0;
00718 }
00719 
00720 
00721 // ------------ method called once each job just after ending the event loop  ------------
00722 void CosmicSplitterValidation::endJob() {
00723         
00724         //std::cout << "totalTracksToAnalyzer: " << totalTracksToAnalyzer_ << std::endl;
00725         std::cout << "golden: " << goldenCtr << ", two tracks: " << twoTracksCtr << ", golden+twotracks: " << goldenPlusTwoTracksCtr << ", tracks+muons cuts: " << _passesTracksPlusMuonsCuts << std::endl;
00726 }
00727 
00728 bool CosmicSplitterValidation::is_gold_muon(const edm::Event& e){
00729         
00730         edm::Handle<reco::MuonCollection> muHandle;
00731         e.getByLabel("STAMuons", muHandle);
00732         const reco::MuonCollection & muons = *(muHandle.product());
00733         // make sure there are 2 muons
00734         if ( 2 != muons.size() ) return false;
00735         
00736         double mudd0=0., mudphi=0., muddsz=0., mudeta=0.;
00737         for ( unsigned int bindex = 0; bindex < muons.size(); ++bindex ) {
00738                 reco::Muon mymuon = muons[bindex];
00739                 // deprecated in 21x (now outerTrack)
00740                 //reco::TrackRef mutrackref = mymuon.standAloneMuon();
00741                 reco::TrackRef mutrackref = mymuon.outerTrack();
00742                 const reco::Track* mutrack = mutrackref.get();
00743                 if (0 == bindex){
00744                         mudd0+=mutrack->d0(); 
00745                         mudphi+=mutrack->phi();
00746                         muddsz+=mutrack->dsz(); 
00747                         mudeta+=mymuon.eta();
00748                 }
00749                 if (1 == bindex){
00750                         mudd0-=mutrack->d0(); 
00751                         mudphi-=mutrack->phi();
00752                         muddsz-=mutrack->dsz(); 
00753                         mudeta-=mymuon.eta();
00754                 }
00755         }
00756         if ((fabs(mudd0)<15.0)&&(fabs(mudphi)<0.045)&&(fabs(muddsz)<20.0)&&(fabs(mudeta)<0.060)) return true;
00757         return false;
00758 }
00759 
00760 
00761 //define this as a plug-in
00762 DEFINE_FWK_MODULE(CosmicSplitterValidation);
/data/refman/pasoursint/CMSSW_6_1_2_SLHC4_patch1/src/Alignment/OfflineValidation/plugins/CosmicSplitterValidation.cc
