diff --git a/DataFormats/MuonReco/interface/ME0Muon.h b/DataFormats/MuonReco/interface/ME0Muon.h index 113e01257d221..50573ef498bcf 100644 --- a/DataFormats/MuonReco/interface/ME0Muon.h +++ b/DataFormats/MuonReco/interface/ME0Muon.h @@ -8,6 +8,9 @@ */ #include "DataFormats/TrackReco/interface/TrackFwd.h" #include "DataFormats/TrackReco/interface/Track.h" +#include "DataFormats/GeometryVector/interface/GlobalPoint.h" +#include "DataFormats/GeometryVector/interface/GlobalVector.h" +#include "DataFormats/Math/interface/AlgebraicROOTObjects.h" #include @@ -53,12 +56,42 @@ namespace reco { double phi() const { return innerTrack_.get()->phi(); } /// pseudorapidity of momentum vector double eta() const { return innerTrack_.get()->eta(); } + + //functions for easy variation of me0muon criteria + /* double xpull() const { return xpull_; } */ + /* double xdiff() const { return xdiff_; } */ + /* double ypull() const { return ypull_; } */ + /* double ydiff() const { return ydiff_; } */ + /* double phidirdiff() const { return phidirdiff_; } */ + + /* void setXpull( const double xpull ) { xpull_ = xpull; } */ + /* void setXdiff( const double xdiff ) { xdiff_ = xdiff; } */ + /* void setYpull( const double ypull ) { ypull_ = ypull; } */ + /* void setYdiff( const double ydiff ) { ydiff_ = ydiff; } */ + /* void setPhidirdiff( const double phidirdiff ) { phidirdiff_ = phidirdiff; } */ + + GlobalPoint globalTrackPosAtSurface() const { return globalTrackPosAtSurface_; } + GlobalVector globalTrackMomAtSurface() const { return globalTrackMomAtSurface_; } + int trackCharge() const { return trackCharge_; } + AlgebraicSymMatrix66 trackCov() const { return trackCov_; } + + void setGlobalTrackPosAtSurface(const GlobalPoint globalTrackPosAtSurface) { globalTrackPosAtSurface_ = globalTrackPosAtSurface; } + void setGlobalTrackMomAtSurface(const GlobalVector globalTrackMomAtSurface) { globalTrackMomAtSurface_ = globalTrackMomAtSurface; } + void setTrackCharge(const int trackCharge) { trackCharge_ = trackCharge; } + void setTrackCov(const AlgebraicSymMatrix66 trackCov) { trackCov_ = trackCov; } private: /// reference to Track reconstructed in the tracker only TrackRef innerTrack_; ME0Segment me0Segment_; int me0segid_; + + GlobalPoint globalTrackPosAtSurface_; + GlobalVector globalTrackMomAtSurface_; + int trackCharge_; + AlgebraicSymMatrix66 trackCov_; + + //double xpull_,ypull_,xdiff_,ydiff_,phidirdiff_; }; } diff --git a/DataFormats/MuonReco/src/classes_def.xml b/DataFormats/MuonReco/src/classes_def.xml index e98c5886dd4bc..cf0f7285c34e1 100755 --- a/DataFormats/MuonReco/src/classes_def.xml +++ b/DataFormats/MuonReco/src/classes_def.xml @@ -184,8 +184,8 @@ - - + + diff --git a/RecoMuon/MuonIdentification/interface/ME0MuonSelector.h b/RecoMuon/MuonIdentification/interface/ME0MuonSelector.h new file mode 100644 index 0000000000000..08e2fcc76de44 --- /dev/null +++ b/RecoMuon/MuonIdentification/interface/ME0MuonSelector.h @@ -0,0 +1,42 @@ +#ifndef RecoMuon_ME0MuonSelectors_h +#define RecoMuon_ME0MuonSelectors_h +// + +#include "DataFormats/MuonReco/interface/ME0Muon.h" + +#include "Geometry/GEMGeometry/interface/ME0Geometry.h" +#include +#include +#include +#include "FWCore/Framework/interface/ESHandle.h" + +#include "TMath.h" +#include + + + + +namespace muon { + /// Selector type + enum SelectionType { + All = 0, // dummy options - always true + VeryLoose = 1, // + Loose = 2, // + Tight = 3, // + }; + + /// a lightweight "map" for selection type string label and enum value + struct SelectionTypeStringToEnum { const char *label; SelectionType value; }; + SelectionType selectionTypeFromString( const std::string &label ); + + /// main GoodMuon wrapper call + bool isGoodMuon( edm::ESHandle me0Geom, const reco::ME0Muon& me0muon, SelectionType type ); + + + /// Specialized isGoodMuon function called from main wrapper + + bool isGoodMuon( edm::ESHandle me0Geom, const reco::ME0Muon& me0muon, double MaxPullX, double MaxDiffX, double MaxPullY, double MaxDiffY, double MaxDiffPhiDir ); + + +} +#endif diff --git a/RecoMuon/MuonIdentification/plugins/ME0MuonSelector.cc b/RecoMuon/MuonIdentification/plugins/ME0MuonSelector.cc new file mode 100644 index 0000000000000..35903870d7e07 --- /dev/null +++ b/RecoMuon/MuonIdentification/plugins/ME0MuonSelector.cc @@ -0,0 +1,121 @@ +/** + *Filter to select me0Muons based on pulls and differences w.r.t. me0Segments + * + * + */ + +#include "RecoMuon/MuonIdentification/interface/ME0MuonSelector.h" + +#include "DataFormats/TrajectoryState/interface/LocalTrajectoryParameters.h" +#include "TrackingTools/AnalyticalJacobians/interface/JacobianCartesianToLocal.h" +#include "TrackingTools/AnalyticalJacobians/interface/JacobianLocalToCartesian.h" + + +namespace muon { +SelectionType selectionTypeFromString( const std::string &label ) +{ + static SelectionTypeStringToEnum selectionTypeStringToEnumMap[] = { + { "All", All }, + { "VeryLoose", VeryLoose }, + { "Loose", Loose }, + { "Tight", Tight }, + { 0, (SelectionType)-1 } + }; + + SelectionType value = (SelectionType)-1; + bool found = false; + for(int i = 0; selectionTypeStringToEnumMap[i].label && (! found); ++i) + if (! strcmp(label.c_str(), selectionTypeStringToEnumMap[i].label)) { + found = true; + value = selectionTypeStringToEnumMap[i].value; + } + + // in case of unrecognized selection type + if (! found) throw cms::Exception("MuonSelectorError") << label << " is not a recognized SelectionType"; + return value; +} +} + + +bool muon::isGoodMuon(edm::ESHandle me0Geom, const reco::ME0Muon& me0muon, SelectionType type) +{ + switch (type) + { + case muon::All: + return true; + break; + case muon::VeryLoose: + return isGoodMuon(me0Geom, me0muon,3,4,20,20,3.14); + break; + case muon::Loose: + return isGoodMuon(me0Geom, me0muon,3,2,3,2,0.5); + break; + case muon::Tight: + return isGoodMuon(me0Geom, me0muon,3,2,3,2,0.15); + break; + default: + return false; + } +} + + + + +bool muon::isGoodMuon(edm::ESHandle me0Geom, const reco::ME0Muon& me0muon, double MaxPullX, double MaxDiffX, double MaxPullY, double MaxDiffY, double MaxDiffPhiDir ) +{ + ME0Segment thisSegment = me0muon.me0segment(); + + ME0DetId id = thisSegment.me0DetId(); + + auto roll = me0Geom->etaPartition(id); + + float zSign = me0muon.globalTrackMomAtSurface().z()/fabs(me0muon.globalTrackMomAtSurface().z()); + if ( zSign * roll->toGlobal(thisSegment.localPosition()).z() < 0 ) return false; + + //GlobalPoint r3FinalReco_glob(r3FinalReco_globv.x(),r3FinalReco_globv.y(),r3FinalReco_globv.z()); + + LocalPoint r3FinalReco = roll->toLocal(me0muon.globalTrackPosAtSurface()); + LocalVector p3FinalReco=roll->toLocal(me0muon.globalTrackMomAtSurface()); + + LocalPoint thisPosition(thisSegment.localPosition()); + LocalVector thisDirection(thisSegment.localDirection().x(),thisSegment.localDirection().y(),thisSegment.localDirection().z()); //FIXME + + //The same goes for the error + AlgebraicMatrix thisCov(4,4,0); + for (int i = 1; i <=4; i++){ + for (int j = 1; j <=4; j++){ + thisCov(i,j) = thisSegment.parametersError()(i,j); + } + } + + ///////////////////////////////////////////////////////////////////////////////////////// + + + LocalTrajectoryParameters ltp(r3FinalReco,p3FinalReco,me0muon.trackCharge()); + JacobianCartesianToLocal jctl(roll->surface(),ltp); + AlgebraicMatrix56 jacobGlbToLoc = jctl.jacobian(); + + AlgebraicMatrix55 Ctmp = (jacobGlbToLoc * me0muon.trackCov()) * ROOT::Math::Transpose(jacobGlbToLoc); + AlgebraicSymMatrix55 C; // I couldn't find any other way, so I resort to the brute force + for(int i=0; i<5; ++i) { + for(int j=0; j<5; ++j) { + C[i][j] = Ctmp[i][j]; + + } + } + + Double_t sigmax = sqrt(C[3][3]+thisSegment.localPositionError().xx() ); + Double_t sigmay = sqrt(C[4][4]+thisSegment.localPositionError().yy() ); + + bool X_MatchFound = false, Y_MatchFound = false, Dir_MatchFound = false; + + + if ( ( (fabs(thisPosition.x()-r3FinalReco.x())/sigmax ) < MaxPullX ) || (fabs(thisPosition.x()-r3FinalReco.x()) < MaxDiffX ) ) X_MatchFound = true; + if ( ( (fabs(thisPosition.y()-r3FinalReco.y())/sigmay ) < MaxPullY ) || (fabs(thisPosition.y()-r3FinalReco.y()) < MaxDiffY ) ) Y_MatchFound = true; + + if ( fabs(me0muon.globalTrackMomAtSurface().phi()-roll->toGlobal(thisSegment.localDirection()).phi()) < MaxDiffPhiDir) Dir_MatchFound = true; + + return (X_MatchFound && Y_MatchFound && Dir_MatchFound); + +} + diff --git a/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.cc b/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.cc index 50367ee747806..0ca01910fc997 100644 --- a/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.cc +++ b/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.cc @@ -39,6 +39,11 @@ ME0SegmentMatcher::ME0SegmentMatcher(const edm::ParameterSet& pas) : iev(0){ produces >(); + X_PULL_CUT = pas.getParameter("maxPullX"); + X_RESIDUAL_CUT = pas.getParameter("maxDiffX"); + Y_PULL_CUT = pas.getParameter("maxPullY"); + Y_RESIDUAL_CUT = pas.getParameter("maxDiffY"); + PHIDIR_RESIDUAL_CUT = pas.getParameter("maxDiffPhiDirection"); } ME0SegmentMatcher::~ME0SegmentMatcher() {} @@ -167,10 +172,16 @@ void ME0SegmentMatcher::produce(edm::Event& ev, const edm::EventSetup& setup) { bool X_MatchFound = false, Y_MatchFound = false, Dir_MatchFound = false; - if ( (fabs(thisPosition.x()-r3FinalReco.x()) < (3.0 * sigmax)) || (fabs(thisPosition.x()-r3FinalReco.x()) < 2.0 ) ) X_MatchFound = true; - if ( (fabs(thisPosition.y()-r3FinalReco.y()) < (3.0 * sigmay)) || (fabs(thisPosition.y()-r3FinalReco.y()) < 2.0 ) ) Y_MatchFound = true; + // if ( (fabs(thisPosition.x()-r3FinalReco.x()) < (3.0 * sigmax)) || (fabs(thisPosition.x()-r3FinalReco.x()) < 2.0 ) ) X_MatchFound = true; + // if ( (fabs(thisPosition.y()-r3FinalReco.y()) < (3.0 * sigmay)) || (fabs(thisPosition.y()-r3FinalReco.y()) < 2.0 ) ) Y_MatchFound = true; - if ( fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi()) < 0.15) Dir_MatchFound = true; + // if ( fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi()) < 0.15) Dir_MatchFound = true; + + + if ( (fabs(thisPosition.x()-r3FinalReco.x()) < (X_PULL_CUT * sigmax)) || (fabs(thisPosition.x()-r3FinalReco.x()) < X_RESIDUAL_CUT ) ) X_MatchFound = true; + if ( (fabs(thisPosition.y()-r3FinalReco.y()) < (Y_PULL_CUT * sigmay)) || (fabs(thisPosition.y()-r3FinalReco.y()) < Y_RESIDUAL_CUT ) ) Y_MatchFound = true; + + if ( fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi()) < PHIDIR_RESIDUAL_CUT) Dir_MatchFound = true; //Check for a Match, and if there is a match, check the delR from the segment, keeping only the closest in MuonCandidate if (X_MatchFound && Y_MatchFound && Dir_MatchFound) { @@ -184,6 +195,20 @@ void ME0SegmentMatcher::produce(edm::Event& ev, const edm::EventSetup& setup) { if (thisDelR < ClosestDelR){ ClosestDelR = thisDelR; MuonCandidate = reco::ME0Muon(thisTrackRef,(*thisSegment),SegmentNumber); + + //Setting the variables for easy me0muon selection later on + // MuonCandidate.setXpull( ( fabs(thisPosition.x()-r3FinalReco.x()) / sigmax) ); + // MuonCandidate.setYpull( ( fabs(thisPosition.y()-r3FinalReco.y()) / sigmay) ); + + // MuonCandidate.setXpull( fabs(thisPosition.x()-r3FinalReco.x()) ); + // MuonCandidate.setYpull( fabs(thisPosition.y()-r3FinalReco.y()) ); + + // MuonCandidate.setPhidirdiff(fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi())); + + MuonCandidate.setGlobalTrackPosAtSurface(r3FinalReco_glob); + MuonCandidate.setGlobalTrackMomAtSurface(p3FinalReco_glob); + MuonCandidate.setTrackCharge(chargeReco); + MuonCandidate.setTrackCov(covFinalReco); } } }//End loop for (auto thisSegment = OurSegments->begin(); thisSegment != OurSegments->end(); ++thisSegment,++SegmentNumber) diff --git a/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.h b/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.h index a4fa10cc8623f..fc09e33a86781 100644 --- a/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.h +++ b/RecoMuon/MuonIdentification/plugins/ME0SegmentMatcher.h @@ -70,7 +70,7 @@ class ME0SegmentMatcher : public edm::EDProducer { int iev; // events through edm::ESHandle me0Geom; - + double X_PULL_CUT, Y_PULL_CUT,X_RESIDUAL_CUT,Y_RESIDUAL_CUT, PHIDIR_RESIDUAL_CUT; }; diff --git a/RecoMuon/MuonIdentification/python/me0SegmentMatcher_cfi.py b/RecoMuon/MuonIdentification/python/me0SegmentMatcher_cfi.py index 932dc48655605..c709ec3357ba1 100644 --- a/RecoMuon/MuonIdentification/python/me0SegmentMatcher_cfi.py +++ b/RecoMuon/MuonIdentification/python/me0SegmentMatcher_cfi.py @@ -1,4 +1,10 @@ import FWCore.ParameterSet.Config as cms -me0SegmentMatching = cms.EDProducer("ME0SegmentMatcher") +me0SegmentMatching = cms.EDProducer("ME0SegmentMatcher", + maxPullX = cms.double (3.0), + maxDiffX = cms.double (4.0), + maxPullY = cms.double (20.0), + maxDiffY = cms.double (20.0), + maxDiffPhiDirection = cms.double (3.14), + ) diff --git a/RecoMuon/MuonIdentification/test/ME0MuonAnalyzer.cc b/RecoMuon/MuonIdentification/test/ME0MuonAnalyzer.cc index 61f790c05fba6..d0c4821457e39 100644 --- a/RecoMuon/MuonIdentification/test/ME0MuonAnalyzer.cc +++ b/RecoMuon/MuonIdentification/test/ME0MuonAnalyzer.cc @@ -31,6 +31,8 @@ #include "TLorentzVector.h" #include "DataFormats/HepMCCandidate/interface/GenParticle.h" +#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h" + //#include "DataFormats/TrackReco/interface/TrackFwd.h" #include "DataFormats/TrackReco/interface/Track.h" @@ -55,10 +57,15 @@ #include "SimTracker/TrackerHitAssociation/interface/TrackerHitAssociator.h" #include "SimTracker/Records/interface/TrackAssociatorRecord.h" +#include "DataFormats/MuonReco/interface/Muon.h" + #include "SimTracker/TrackAssociation/plugins/ParametersDefinerForTPESProducer.h" #include "SimTracker/TrackAssociation/plugins/CosmicParametersDefinerForTPESProducer.h" #include "CommonTools/CandAlgos/interface/GenParticleCustomSelector.h" +//#include "CommonTools/CandAlgos/interface/TrackingParticleCustomSelector.h" + +#include "RecoMuon/MuonIdentification/plugins/ME0MuonSelector.cc" #include "Fit/FitResult.h" #include "TF1.h" @@ -93,6 +100,7 @@ #include #include +#include class ME0MuonAnalyzer : public edm::EDAnalyzer { public: @@ -111,9 +119,12 @@ class ME0MuonAnalyzer : public edm::EDAnalyzer { virtual void analyze(const edm::Event&, const edm::EventSetup&); - virtual void endJob(); + //virtual void endJob(); //virtual void beginJob(const edm::EventSetup&); - void beginJob(); + //void beginJob(); + void beginRun(edm::Run const&, edm::EventSetup const&); + void endRun(edm::Run const&, edm::EventSetup const&); + //For Track Association @@ -126,45 +137,85 @@ class ME0MuonAnalyzer : public edm::EDAnalyzer { //edm::InputTag associatormap; bool UseAssociators; bool RejectEndcapMuons; - const TrackAssociatorByChi2* associatorByChi2; + //const TrackAssociatorByChi2* associatorByChi2; + //const TrackAssociatorByHits* associatorByHits; std::vector associators; std::vector associator; std::vector label; - GenParticleCustomSelector gpSelector; + //GenParticleCustomSelector gpSelector; + //TrackingParticleCustomSelector gpSelector; //std::string parametersDefiner; TString histoFolder; + TString me0MuonSelector; TFile* histoFile; TH1F *Candidate_Eta; TH1F *Mass_h; TH1F *Segment_Eta; TH1F *Segment_Phi; TH1F *Segment_R; TH2F *Segment_Pos; TH1F *Rechit_Eta; TH1F *Rechit_Phi; TH1F *Rechit_R; TH2F *Rechit_Pos; TH1F *GenMuon_Phi; TH1F *GenMuon_R; TH2F *GenMuon_Pos; TH1F *Track_Eta; TH1F *Track_Pt; TH1F *ME0Muon_Eta; TH1F *ME0Muon_Pt; TH1F *CheckME0Muon_Eta; - TH1F *ME0Muon_Cuts_Eta_5_10; TH1F *ME0Muon_Cuts_Eta_10_20; TH1F *ME0Muon_Cuts_Eta_20_40; TH1F *ME0Muon_Cuts_Eta_40; - TH1F *UnmatchedME0Muon_Eta; TH1F *UnmatchedME0Muon_Pt; TH1F *Chi2UnmatchedME0Muon_Eta; - TH1F *UnmatchedME0Muon_Cuts_Eta_5_10; TH1F *UnmatchedME0Muon_Cuts_Eta_10_20; TH1F *UnmatchedME0Muon_Cuts_Eta_20_40; TH1F *UnmatchedME0Muon_Cuts_Eta_40; + TH1F *ME0Muon_SmallBins_Pt; TH1F *ME0Muon_VariableBins_Pt; + TH1F *ME0Muon_Cuts_Eta_5_10; TH1F *ME0Muon_Cuts_Eta_9_11; TH1F *ME0Muon_Cuts_Eta_10_50; TH1F *ME0Muon_Cuts_Eta_50_100; TH1F *ME0Muon_Cuts_Eta_100; + TH1F *UnmatchedME0Muon_Eta; TH1F *UnmatchedME0Muon_Pt; TH1F *UnmatchedME0Muon_Window_Pt; TH1F *Chi2UnmatchedME0Muon_Eta; + TH1F *Chi2UnmatchedME0Muon_Pt; TH1F *Chi2UnmatchedME0Muon_SmallBins_Pt; TH1F *Chi2UnmatchedME0Muon_VariableBins_Pt; + TH1F *UnmatchedME0Muon_SmallBins_Pt; TH1F *UnmatchedME0Muon_VariableBins_Pt; + TH1F *UnmatchedME0Muon_Cuts_Eta_5_10; TH1F *UnmatchedME0Muon_Cuts_Eta_9_11; TH1F *UnmatchedME0Muon_Cuts_Eta_10_50; TH1F *UnmatchedME0Muon_Cuts_Eta_50_100; TH1F *UnmatchedME0Muon_Cuts_Eta_100; TH1F *TracksPerSegment_h; TH2F *TracksPerSegment_s; TProfile *TracksPerSegment_p; TH2F *ClosestDelR_s; TProfile *ClosestDelR_p; TH2F *PtDiff_s; TProfile *PtDiff_p; TH1F *PtDiff_h; TH1F *QOverPtDiff_h; TH1F *PtDiff_rms; TH1F *PtDiff_gaus_narrow; TH1F *PtDiff_gaus_wide; + TH2F *StandalonePtDiff_s; TProfile *StandalonePtDiff_p; TH1F *StandalonePtDiff_h; TH1F *StandaloneQOverPtDiff_h; TH1F *StandalonePtDiff_rms; TH1F *StandalonePtDiff_gaus_narrow; TH1F *StandalonePtDiff_gaus_wide; TH1F *PtDiff_gaus_5_10; TH1F *PtDiff_gaus_10_20; TH1F *PtDiff_gaus_20_40; TH1F *PtDiff_gaus_40; + TH1F *StandalonePtDiff_gaus; TH1F *VertexDiff_h; TH2F *PDiff_s; TProfile *PDiff_p; TH1F *PDiff_h; TH2F *PtDiff_s_5_10; TH2F *PtDiff_s_10_20; TH2F *PtDiff_s_20_40; TH2F *PtDiff_s_40; TH1F *FakeTracksPerSegment_h; TH2F *FakeTracksPerSegment_s; TProfile *FakeTracksPerSegment_p; TH1F *FakeTracksPerAssociatedSegment_h; TH2F *FakeTracksPerAssociatedSegment_s; TProfile *FakeTracksPerAssociatedSegment_p; TH1F *GenMuon_Eta; TH1F *GenMuon_Pt; TH1F *MatchedME0Muon_Eta; TH1F *MatchedME0Muon_Pt; TH1F *Chi2MatchedME0Muon_Eta; TH1F *Chi2MatchedME0Muon_Pt; - TH1F *MatchedME0Muon_Eta_5_10; TH1F *MatchedME0Muon_Eta_10_20; TH1F *MatchedME0Muon_Eta_20_40; TH1F *MatchedME0Muon_Eta_40; - TH1F *GenMuon_Eta_5_10; TH1F *GenMuon_Eta_10_20; TH1F *GenMuon_Eta_20_40; TH1F *GenMuon_Eta_40; + TH1F *GenMuon_SmallBins_Pt; TH1F *MatchedME0Muon_SmallBins_Pt; TH1F *Chi2MatchedME0Muons_Pt; TH1F *Chi2MatchedME0Muon_SmallBins_Pt; + TH1F *GenMuon_VariableBins_Pt; TH1F *MatchedME0Muon_VariableBins_Pt; TH1F *Chi2MatchedME0Muon_VariableBins_Pt; + TH1F *TPMuon_Eta; TH1F *TPMuon_SmallBins_Pt; TH1F *TPMuon_Pt; TH1F *TPMuon_VariableBins_Pt; + TH1F *MatchedME0Muon_Eta_5_10; TH1F *MatchedME0Muon_Eta_9_11; TH1F *MatchedME0Muon_Eta_10_20; TH1F *MatchedME0Muon_Eta_20_40; TH1F *MatchedME0Muon_Eta_40; + TH1F *Chi2MatchedME0Muon_Eta_5_10; TH1F *Chi2MatchedME0Muon_Eta_9_11; TH1F *Chi2MatchedME0Muon_Eta_10_20; TH1F *Chi2MatchedME0Muon_Eta_20_40; TH1F *Chi2MatchedME0Muon_Eta_40; + TH1F *TPMuon_Eta_5_10; TH1F *TPMuon_Eta_10_20; TH1F *TPMuon_Eta_20_40; TH1F *TPMuon_Eta_40; + TH1F *GenMuon_Eta_5_10; TH1F *GenMuon_Eta_9_11; TH1F *GenMuon_Eta_10_20; TH1F *GenMuon_Eta_20_40; TH1F *GenMuon_Eta_40; TH1F *MuonRecoEff_Eta; TH1F *MuonRecoEff_Pt; TH1F *Chi2MuonRecoEff_Eta; - TH1F *MuonRecoEff_Eta_5_10; TH1F *MuonRecoEff_Eta_10_20; TH1F *MuonRecoEff_Eta_20_40; TH1F *MuonRecoEff_Eta_40; + TH1F *MuonRecoEff_Eta_5_10; TH1F *MuonRecoEff_Eta_9_11; TH1F *MuonRecoEff_Eta_10_20; TH1F *MuonRecoEff_Eta_20_40; TH1F *MuonRecoEff_Eta_40; + TH1F *Chi2MuonRecoEff_Eta_5_10; TH1F *Chi2MuonRecoEff_Eta_9_11; TH1F *Chi2MuonRecoEff_Eta_10_20; TH1F *Chi2MuonRecoEff_Eta_20_40; TH1F *Chi2MuonRecoEff_Eta_40; TH1F *FakeRate_Eta; TH1F *FakeRate_Pt; TH1F *FakeRate_Eta_PerEvent; TH1F *Chi2FakeRate_Eta; - TH1F *FakeRate_Eta_5_10; TH1F *FakeRate_Eta_10_20; TH1F *FakeRate_Eta_20_40; TH1F *FakeRate_Eta_40; + + TH1F *Chi2FakeRate_WideBinning_Eta; + TH1F *Chi2FakeRate_WidestBinning_Eta; + TH1F *FakeRate_WideBinning_Eta; + TH1F *FakeRate_WidestBinning_Eta; + TH1F *UnmatchedME0Muon_Cuts_WideBinning_Eta; + TH1F *UnmatchedME0Muon_Cuts_WidestBinning_Eta; + TH1F *ME0Muon_Cuts_WideBinning_Eta; + TH1F *ME0Muon_Cuts_WidestBinning_Eta; + TH1F *Chi2UnmatchedME0Muon_WideBinning_Eta; + TH1F *Chi2UnmatchedME0Muon_WidestBinning_Eta; + TH1F *TPMuon_WideBinning_Eta; + TH1F *TPMuon_WidestBinning_Eta; + TH1F *GenMuon_WideBinning_Eta; + TH1F *GenMuon_WidestBinning_Eta; + TH1F *MatchedME0Muon_WideBinning_Eta; + TH1F *MatchedME0Muon_WidestBinning_Eta; + TH1F *Chi2MatchedME0Muon_WideBinning_Eta; + TH1F *Chi2MatchedME0Muon_WidestBinning_Eta; + TH1F *MuonRecoEff_WideBinning_Eta; + TH1F *MuonRecoEff_WidestBinning_Eta; + TH1F *Chi2MuonRecoEff_WideBinning_Eta; + TH1F *Chi2MuonRecoEff_WidestBinning_Eta; + + + TH1F *FakeRate_Eta_5_10; TH1F *FakeRate_Eta_9_11; TH1F *FakeRate_Eta_10_50; TH1F *FakeRate_Eta_50_100; TH1F *FakeRate_Eta_100; TH1F *MuonAllTracksEff_Eta; TH1F *MuonAllTracksEff_Pt; TH1F *MuonUnmatchedTracksEff_Eta; TH1F *MuonUnmatchedTracksEff_Pt; TH1F *FractionMatched_Eta; + TH1F *StandaloneMuonRecoEff_Eta; TH1F *StandaloneMuonRecoEff_WideBinning_Eta; TH1F *StandaloneMuonRecoEff_WidestBinning_Eta; TH1F *UnmatchedME0Muon_Cuts_Eta;TH1F *ME0Muon_Cuts_Eta; - + TH1F *StandaloneMatchedME0Muon_Eta; TH1F *StandaloneMatchedME0Muon_WideBinning_Eta; TH1F *StandaloneMatchedME0Muon_WidestBinning_Eta; TH1F *DelR_Segment_GenMuon; TH1F *SegPosDirPhiDiff_True_h; TH1F *SegPosDirEtaDiff_True_h; TH1F *SegPosDirPhiDiff_All_h; TH1F *SegPosDirEtaDiff_All_h; @@ -175,21 +226,35 @@ class ME0MuonAnalyzer : public edm::EDAnalyzer { TH1F *XDiff_h; TH1F *YDiff_h; TH1F *XPull_h; TH1F *YPull_h; + TH1F *DelR_Window_Under5; TH1F *Pt_Window_Under5; + TH1F *DelR_Track_Window_Under5; TH1F *Pt_Track_Window_Under5; TH1F *Pt_Track_Window; + TH1F *DelR_Track_Window_Failed_Under5; TH1F *Pt_Track_Window_Failed_Under5; TH1F *Pt_Track_Window_Failed; + + TH1F *FailedTrack_Window_XPull; TH1F *FailedTrack_Window_YPull; TH1F *FailedTrack_Window_PhiDiff; + TH1F *FailedTrack_Window_XDiff; TH1F *FailedTrack_Window_YDiff; TH1F *NormChi2_h; TH1F *NormChi2Prob_h; TH2F *NormChi2VsHits_h; TH2F *chi2_vs_eta_h; TH1F *AssociatedChi2_h; TH1F *AssociatedChi2_Prob_h; + TH1F *PreMatch_TP_R; TH1F *PostMatch_TP_R; TH1F *PostMatch_BX0_TP_R; + + TH2F *UnmatchedME0Muon_ScatterPlot; + double FakeRatePtCut, MatchingWindowDelR; double Nevents; + TH1F *Nevents_h; + //Removing this }; ME0MuonAnalyzer::ME0MuonAnalyzer(const edm::ParameterSet& iConfig) { + std::cout<<"Contructor"<("HistoFile").c_str(), "recreate"); histoFolder = iConfig.getParameter("HistoFolder").c_str(); + me0MuonSelector = iConfig.getParameter("ME0MuonSelectionType").c_str(); RejectEndcapMuons = iConfig.getParameter< bool >("RejectEndcapMuons"); UseAssociators = iConfig.getParameter< bool >("UseAssociators"); @@ -201,114 +266,187 @@ ME0MuonAnalyzer::ME0MuonAnalyzer(const edm::ParameterSet& iConfig) UseAssociators = iConfig.getParameter< bool >("UseAssociators"); associators = iConfig.getParameter< std::vector >("associators"); - label = iConfig.getParameter< std::vector >("label"), + label = iConfig.getParameter< std::vector >("label"); + + // gpSelector = GenParticleCustomSelector(iConfig.getParameter("ptMinGP"), + // iConfig.getParameter("minRapidityGP"), + // iConfig.getParameter("maxRapidityGP"), + // iConfig.getParameter("tipGP"), + // iConfig.getParameter("lipGP"), + // iConfig.getParameter("chargedOnlyGP"), + // iConfig.getParameter("statusGP"), + // iConfig.getParameter >("pdgIdGP")); + + // gpSelector = TrackingParticleCustomSelector(iConfig.getParameter("ptMinGP"), + // iConfig.getParameter("minRapidityGP"), + // iConfig.getParameter("maxRapidityGP"), + // iConfig.getParameter("tipGP"), + // iConfig.getParameter("lipGP"), + // iConfig.getParameter("chargedOnlyGP"), + // iConfig.getParameter("statusGP"), + // iConfig.getParameter >("pdgIdGP")); + // // + //parametersDefiner =iConfig.getParameter("parametersDefiner"); + + std::cout<<"Contructor end"<("ptMinGP"), - iConfig.getParameter("minRapidityGP"), - iConfig.getParameter("maxRapidityGP"), - iConfig.getParameter("tipGP"), - iConfig.getParameter("lipGP"), - iConfig.getParameter("chargedOnlyGP"), - iConfig.getParameter("statusGP"), - iConfig.getParameter >("pdgIdGP")); - //parametersDefiner =iConfig.getParameter("parametersDefiner"); -} +//void ME0MuonAnalyzer::beginJob(const edm::EventSetup& iSetup) +void ME0MuonAnalyzer::beginRun(edm::Run const&, edm::EventSetup const& iSetup) { +//void ME0MuonAnalyzer::beginJob() +//{ + std::cout<<"At start of begin run"< theAssociator; + std::cout<<"HERE NOW"<().get(associators[w],theAssociator); + associator.push_back( theAssociator.product() ); + } + } } @@ -406,8 +611,8 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup { - //std::cout<<"ANALYZER"<Fill(1); using namespace edm; //run_ = (int)iEvent.id().run(); @@ -423,7 +628,8 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup // iEvent.getByLabel ("me0SegmentMatcher", OurMuons); Handle > OurCandidates; - iEvent.getByLabel > ("me0MuonConverter", OurCandidates); + //iEvent.getByLabel > ("me0MuonConverter", OurCandidates); + iEvent.getByLabel > ("me0MuonConverting", OurCandidates); //Handle > OurSegments; //iEvent.getByLabel >("me0SegmentProducer", OurSegments); @@ -435,13 +641,17 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup unsigned int gensize=genParticles->size(); + Handle trackingParticles; + iEvent.getByLabel("mix","MergedTrackTruth", trackingParticles); + + if (RejectEndcapMuons){ //Section to turn off signal muons in the endcaps, to approximate a nu gun for(unsigned int i=0; i 1.9 ) { - std::cout<<"Found a signal muon outside the barrel, exiting the function"< generalTracks; iEvent.getByLabel ("generalTracks", generalTracks); - + //std::cout<<"About to get me0muons"< > OurMuons; - iEvent.getByLabel > ("me0SegmentMatcher", OurMuons); + iEvent.getByLabel > ("me0SegmentMatching", OurMuons); Handle OurSegments; iEvent.getByLabel("me0Segments","",OurSegments); @@ -475,20 +685,6 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup //For Track Association: - //std::cout<<"ON BEGIN JOB:"< theAssociator; - //std::cout<<"associators size = "<().get(associators[w],theAssociator); - //std::cout<<"On step "< parametersDefinerTP; @@ -500,8 +696,10 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup std::vector IsMatched; std::vector SegIdForMatch; + //std::cout<<"Looping on me0Muons"<::const_iterator thisMuon = OurMuons->begin(); thisMuon != OurMuons->end(); ++thisMuon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; IsMatched.push_back(false); SegIdForMatch.push_back(-1); } @@ -688,30 +886,492 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup // // } // } - //Track Association by Chi2: + + std::vector MatchedSegIds; + + + for(unsigned int i=0; i 2.0) && (fabs(CurrentParticle.eta()) < 2.8) ) GenMuon_Eta->Fill(CurrentParticle.eta()); + ////std::cout<<"Mother's ID is: "< ReferenceTrackPt; + + double VertexDiff=-1,PtDiff=-1,QOverPtDiff=-1,PDiff=-1; + + + //std::cout<<"Size = "<size()<::const_iterator thisMuon = OurMuons->begin(); + thisMuon != OurMuons->end(); ++thisMuon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; + TrackRef tkRef = thisMuon->innerTrack(); + SegIdForMatch.push_back(thisMuon->me0segid()); + thisDelR = reco::deltaR(CurrentParticle,*tkRef); + ReferenceTrackPt.push_back(tkRef->pt()); + + + if (( tkRef->pt() > FakeRatePtCut ) ){ + if (thisDelR < MatchingWindowDelR ){ + if (tkRef->pt() < 5.0){ + DelR_Window_Under5->Fill(thisDelR); + Pt_Window_Under5->Fill(tkRef->pt()); + } + if (thisDelR < LowestDelR){ + LowestDelR = thisDelR; + //if (fabs(tkRef->pt() - CurrentParticle.pt())/CurrentParticle.pt() < 0.50) MatchedID = ME0MuonID; + MatchedID = ME0MuonID; + VertexDiff = fabs(tkRef->vz()-CurrentParticle.vz()); + QOverPtDiff = ( (tkRef->charge() /tkRef->pt()) - (CurrentParticle.charge()/CurrentParticle.pt() ) )/ (CurrentParticle.charge()/CurrentParticle.pt() ); + PtDiff = (tkRef->pt() - CurrentParticle.pt())/CurrentParticle.pt(); + PDiff = (tkRef->p() - CurrentParticle.p())/CurrentParticle.p(); + } + } + } + + ME0MuonID++; + + } + + for (std::vector::const_iterator thisTrack = generalTracks->begin(); + thisTrack != generalTracks->end();++thisTrack){ + //TrackRef tkRef = thisTrack->innerTrack(); + thisDelR = reco::deltaR(CurrentParticle,*thisTrack); + + if ((thisTrack->pt() > FakeRatePtCut ) ){ + if (thisDelR < MatchingWindowDelR ){ + if (thisTrack->pt() < 5.0){ + DelR_Track_Window_Under5->Fill(thisDelR); + Pt_Track_Window_Under5->Fill(thisTrack->pt()); + } + Pt_Track_Window->Fill(thisTrack->pt()); + } + } + } + if (MatchedID == -1){ + + for (std::vector::const_iterator thisTrack = generalTracks->begin(); + thisTrack != generalTracks->end();++thisTrack){ + //TrackRef tkRef = thisTrack->innerTrack(); + thisDelR = reco::deltaR(CurrentParticle,*thisTrack); + + + if ( (thisTrack->pt() > FakeRatePtCut ) && (TMath::Abs(thisTrack->eta()) < 2.8) && (TMath::Abs(thisTrack->eta()) > 2.0) ) { + if (thisDelR < MatchingWindowDelR ){ + if (thisTrack->pt() < 5.0){ + DelR_Track_Window_Failed_Under5->Fill(thisDelR); + Pt_Track_Window_Failed_Under5->Fill(thisTrack->pt()); + } + // if ( (thisTrack->pt() > 9.5) && (thisTrack->pt() < 10.5) ){ + // //if ( (thisTrack->pt() > 0) ){ + // ////////////////Start propagation from ME0SegmentMatcher + // float zSign = thisTrack->pz()/fabs(thisTrack->pz()); + + // float zValue = 526.75 * zSign; + + // Plane *plane = new Plane(Surface::PositionType(0,0,zValue),Surface::RotationType()); + + // //Getting the initial variables for propagation + + // int chargeReco = thisTrack->charge(); + // GlobalVector p3reco, r3reco; + + // p3reco = GlobalVector(thisTrack->outerPx(), thisTrack->outerPy(), thisTrack->outerPz()); + // r3reco = GlobalVector(thisTrack->outerX(), thisTrack->outerY(), thisTrack->outerZ()); + + // AlgebraicSymMatrix66 covReco; + // //This is to fill the cov matrix correctly + // AlgebraicSymMatrix55 covReco_curv; + // covReco_curv = thisTrack->outerStateCovariance(); + // FreeTrajectoryState initrecostate = getFTS(p3reco, r3reco, chargeReco, covReco_curv, &*bField); + // getFromFTS(initrecostate, p3reco, r3reco, chargeReco, covReco); + + // //Now we propagate and get the propagated variables from the propagated state + // SteppingHelixStateInfo startrecostate(initrecostate); + // SteppingHelixStateInfo lastrecostate; + + // const SteppingHelixPropagator* ThisshProp = + // dynamic_cast(&*shProp); + + // lastrecostate = ThisshProp->propagate(startrecostate, *plane); + + // FreeTrajectoryState finalrecostate; + // lastrecostate.getFreeState(finalrecostate); + + // AlgebraicSymMatrix66 covFinalReco; + // GlobalVector p3FinalReco_glob, r3FinalReco_globv; + // getFromFTS(finalrecostate, p3FinalReco_glob, r3FinalReco_globv, chargeReco, covFinalReco); + + + // //To transform the global propagated track to local coordinates + // int SegmentNumber = 0; + + // reco::ME0Muon MuonCandidate; + + // for (auto thisSegment = OurSegments->begin(); thisSegment != OurSegments->end(); + // ++thisSegment,++SegmentNumber){ + // ME0DetId id = thisSegment->me0DetId(); + + // auto roll = me0Geom->etaPartition(id); + + // if ( zSign * roll->toGlobal(thisSegment->localPosition()).z() < 0 ) continue; + + // GlobalPoint r3FinalReco_glob(r3FinalReco_globv.x(),r3FinalReco_globv.y(),r3FinalReco_globv.z()); + + // LocalPoint r3FinalReco = roll->toLocal(r3FinalReco_glob); + // LocalVector p3FinalReco=roll->toLocal(p3FinalReco_glob); + + // LocalPoint thisPosition(thisSegment->localPosition()); + // LocalVector thisDirection(thisSegment->localDirection().x(),thisSegment->localDirection().y(),thisSegment->localDirection().z()); //FIXME + + // //The same goes for the error + // AlgebraicMatrix thisCov(4,4,0); + // for (int i = 1; i <=4; i++){ + // for (int j = 1; j <=4; j++){ + // thisCov(i,j) = thisSegment->parametersError()(i,j); + // } + // } + + // ///////////////////////////////////////////////////////////////////////////////////////// + + + // LocalTrajectoryParameters ltp(r3FinalReco,p3FinalReco,chargeReco); + // JacobianCartesianToLocal jctl(roll->surface(),ltp); + // AlgebraicMatrix56 jacobGlbToLoc = jctl.jacobian(); + + // AlgebraicMatrix55 Ctmp = (jacobGlbToLoc * covFinalReco) * ROOT::Math::Transpose(jacobGlbToLoc); + // AlgebraicSymMatrix55 C; // I couldn't find any other way, so I resort to the brute force + // for(int i=0; i<5; ++i) { + // for(int j=0; j<5; ++j) { + // C[i][j] = Ctmp[i][j]; + + // } + // } + + // Double_t sigmax = sqrt(C[3][3]+thisSegment->localPositionError().xx() ); + // Double_t sigmay = sqrt(C[4][4]+thisSegment->localPositionError().yy() ); + + + // if( reco::deltaR(roll->toGlobal(thisSegment->localPosition()),*thisTrack) < 0.3){ + + // //if( reco::deltaR(roll->toGlobal(thisSegment->localPosition()),r3FinalReco) < 0.3){ + + // FailedTrack_Window_XPull->Fill(fabs(thisPosition.x()-r3FinalReco.x()) / (1.0 * sigmax)); + // FailedTrack_Window_YPull->Fill(fabs(thisPosition.y()-r3FinalReco.y()) / (1.0 * sigmay)); + + // FailedTrack_Window_XDiff->Fill(fabs(thisPosition.x()-r3FinalReco.x())); + // FailedTrack_Window_YDiff->Fill(fabs(thisPosition.y()-r3FinalReco.y())); + + + // //Calculating if there's a match in distance, if so, we also plot the direction difference + // //bool X_MatchFound = false, Y_MatchFound = false; + // //if ( (fabs(thisPosition.x()-r3FinalReco.x()) < (3.0 * sigmax)) || (fabs(thisPosition.x()-r3FinalReco.x()) < 2.0 ) ) X_MatchFound = true; + // //if ( (fabs(thisPosition.y()-r3FinalReco.y()) < (3.0 * sigmay)) || (fabs(thisPosition.y()-r3FinalReco.y()) < 2.0 ) ) Y_MatchFound = true; + + // //if (X_MatchFound && Y_MatchFound) FailedTrack_Window_PhiDiff->Fill(fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi())); + // FailedTrack_Window_PhiDiff->Fill(fabs(p3FinalReco_glob.phi()-roll->toGlobal(thisSegment->localDirection()).phi())); + + // ////////////////End propagation from ME0SegmentMatcher + // } + + // } + // } + Pt_Track_Window_Failed->Fill(thisTrack->pt()); + } + } + } + } + + if (MatchedID != -1){ + IsMatched[MatchedID] = true; + + if ((CurrentParticle.pt() >FakeRatePtCut) ){ + MatchedME0Muon_Eta->Fill(fabs(CurrentParticle.eta())); + if ( (TMath::Abs(CurrentParticle.eta()) > 2.0) && (TMath::Abs(CurrentParticle.eta()) < 2.8) ) { + MatchedME0Muon_Pt->Fill(CurrentParticle.pt()); + MatchedME0Muon_SmallBins_Pt->Fill(CurrentParticle.pt()); + MatchedME0Muon_VariableBins_Pt->Fill(CurrentParticle.pt()); + } + + + MatchedME0Muon_WideBinning_Eta->Fill(fabs(CurrentParticle.eta())); + MatchedME0Muon_WidestBinning_Eta->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) MatchedME0Muon_Eta_5_10->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 9.0) && (CurrentParticle.pt() <= 11.0) ) MatchedME0Muon_Eta_9_11->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) MatchedME0Muon_Eta_10_20->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) MatchedME0Muon_Eta_20_40->Fill(fabs(CurrentParticle.eta())); + if ( CurrentParticle.pt() > 40.0) MatchedME0Muon_Eta_40->Fill(fabs(CurrentParticle.eta())); + + + + VertexDiff_h->Fill(VertexDiff); + PtDiff_h->Fill(PtDiff); + QOverPtDiff_h->Fill(QOverPtDiff); + //PtDiff_s->Fill(CurrentParticle.eta(),PtDiff); + PtDiff_s->Fill(CurrentParticle.eta(),QOverPtDiff); + if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) PtDiff_s_5_10->Fill(CurrentParticle.eta(),QOverPtDiff); + if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) PtDiff_s_10_20->Fill(CurrentParticle.eta(),QOverPtDiff); + if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) PtDiff_s_20_40->Fill(CurrentParticle.eta(),QOverPtDiff); + if ( CurrentParticle.pt() > 40.0) PtDiff_s_40->Fill(CurrentParticle.eta(),QOverPtDiff); + PtDiff_p->Fill(CurrentParticle.eta(),PtDiff); + + PDiff_h->Fill(PDiff); + PDiff_s->Fill(CurrentParticle.eta(),PDiff); + PDiff_p->Fill(CurrentParticle.eta(),PDiff); + } + MatchedSegIds.push_back(SegIdForMatch[MatchedID]); + + // if ( ((CurrentParticle.eta()) > 2.4) && ((CurrentParticle.eta()) < 3.8) ) { + // //MatchedME0Muon_Pt->Fill(CurrentParticle.pt()); + // //MatchedME0Muon_Pt->Fill(.pt()); + + // } + } + + + if ( (CurrentParticle.pt() >FakeRatePtCut) ){ + GenMuon_Eta->Fill(fabs(CurrentParticle.eta())); + GenMuon_WideBinning_Eta->Fill(fabs(CurrentParticle.eta())); + GenMuon_WidestBinning_Eta->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) GenMuon_Eta_5_10->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 9.0) && (CurrentParticle.pt() <= 11.0) ) GenMuon_Eta_9_11->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) GenMuon_Eta_10_20->Fill(fabs(CurrentParticle.eta())); + if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) GenMuon_Eta_20_40->Fill(fabs(CurrentParticle.eta())); + if ( CurrentParticle.pt() > 40.0) GenMuon_Eta_40->Fill(fabs(CurrentParticle.eta())); + GenMuon_Phi->Fill(CurrentParticle.phi()); + if ( (fabs(CurrentParticle.eta()) > 2.0) && (fabs(CurrentParticle.eta()) < 2.8) ) { + GenMuon_SmallBins_Pt->Fill(CurrentParticle.pt()); + GenMuon_VariableBins_Pt->Fill(CurrentParticle.pt()); + GenMuon_Pt->Fill(CurrentParticle.pt()); + } + } + + } + } + + +//Diff studies for gen matching + + // std::cout<<"Doing first propagation"< 2.8 )) continue; + + // float zSign = CurrentParticle.pz()/fabs(CurrentParticle.pz()); + + // float zValue = 526.75 * zSign; + // Plane *plane = new Plane(Surface::PositionType(0,0,zValue),Surface::RotationType()); + // TLorentzVector Momentum; + // Momentum.SetPtEtaPhiM(CurrentParticle.pt() + // ,CurrentParticle.eta() + // ,CurrentParticle.phi() + // ,CurrentParticle.mass()); + // GlobalVector p3gen(Momentum.Px(), Momentum.Py(), Momentum.Pz()); + // GlobalVector r3gen = GlobalVector(CurrentParticle.vertex().x() + // ,CurrentParticle.vertex().y() + // ,CurrentParticle.vertex().z()); + + // AlgebraicSymMatrix66 covGen = AlgebraicMatrixID(); + // covGen *= 1e-20; // initialize to sigma=1e-10 .. should get overwhelmed by MULS + // AlgebraicSymMatrix66 covFinal; + // int chargeGen = CurrentParticle.charge(); + + // //Propagation + // FreeTrajectoryState initstate = getFTS(p3gen, r3gen, chargeGen, covGen, &*bField); + + // SteppingHelixStateInfo startstate(initstate); + // SteppingHelixStateInfo laststate; + + // const SteppingHelixPropagator* ThisshProp = + // dynamic_cast(&*shProp); + + // laststate = ThisshProp->propagate(startstate, *plane); + + // FreeTrajectoryState finalstate; + // laststate.getFreeState(finalstate); + + // GlobalVector p3Final, r3Final; + // getFromFTS(finalstate, p3Final, r3Final, chargeGen, covFinal); + + + // int ME0MuonID = 0; + + // for (std::vector::const_iterator thisMuon = OurMuons->begin(); + // thisMuon != OurMuons->end(); ++thisMuon){ + + // TrackRef tkRef = thisMuon->innerTrack(); + // SegIdForMatch.push_back(thisMuon->me0segid()); + + // ME0Segment Seg = thisMuon->me0segment(); + // ME0DetId id =Seg.me0DetId(); + // auto roll = me0Geom->etaPartition(id); + + // double DirectionPull, DirectionPullNum, DirectionPullDenom; + + // //Computing the sigma for the track direction + // Double_t mag_track = p3Final.perp(); + // //Double_t phi_track = p3Final.phi(); + + // //Double_t dmagdx_track = p3Final.x()/mag_track; + // //Double_t dmagdy_track = p3Final.y()/mag_track; + // Double_t dphidx_track = -p3Final.y()/(mag_track*mag_track); + // Double_t dphidy_track = p3Final.x()/(mag_track*mag_track); + // Double_t sigmaphi_track = sqrt( dphidx_track*dphidx_track*covFinal(3,3)+ + // dphidy_track*dphidy_track*covFinal(4,4)+ + // dphidx_track*dphidy_track*2*covFinal(3,4) ); + + // DirectionPullNum = p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi(); + // DirectionPullDenom = sqrt( pow(roll->toGlobal(Seg.localPosition()).phi(),2) + pow(sigmaphi_track,2) ); + // DirectionPull = DirectionPullNum / DirectionPullDenom; + + + // if (IsMatched[ME0MuonID]){ + // SegGenDirPhiDiff_True_h->Fill(p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegGenDirEtaDiff_True_h->Fill(p3Final.eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // SegGenDirPhiPull_True_h->Fill(DirectionPull); + // } + + // if ((zSign * roll->toGlobal(Seg.localDirection()).z()) > 0 ){ + // SegGenDirPhiDiff_All_h->Fill(p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegGenDirPhiPull_All_h->Fill(DirectionPull); + + // SegGenDirEtaDiff_All_h->Fill(p3Final.eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // } + // ME0MuonID++; + // } + // } + // } + + //Del R study =========================== + for(unsigned int i=0; i::const_iterator thisMuon = OurMuons->begin(); + thisMuon != OurMuons->end(); ++thisMuon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; + TrackRef tkRef = thisMuon->innerTrack(); + thisDelR = reco::deltaR(CurrentParticle,*tkRef); + if (thisDelR < LowestDelR) LowestDelR = thisDelR; + } + + ClosestDelR_s->Fill(CurrentParticle.eta(), LowestDelR); + ClosestDelR_p->Fill(CurrentParticle.eta(), LowestDelR); + } + } + + //==================================== + + // -----Finally, we loop over all the ME0Muons in the event + // -----Before, we plotted the gen muon pt and eta for the efficiency plot of matches + // -----Now, each time a match failed, we plot the ME0Muon pt and eta + int ME0MuonID = 0; + for (std::vector::const_iterator thisMuon = OurMuons->begin(); + thisMuon != OurMuons->end(); ++thisMuon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; + TrackRef tkRef = thisMuon->innerTrack(); + + + + + //Moved resolution stuff here, only calculate resolutions for matched muons! + + + if (!IsMatched[ME0MuonID]){ + + UnmatchedME0Muon_Eta->Fill(fabs(tkRef->eta())); + + + + if ((tkRef->pt() > FakeRatePtCut) && (TMath::Abs(tkRef->eta()) < 2.8) ) { + if ( (tkRef->pt() > 5.0) && (tkRef->pt() <= 10.0) ) UnmatchedME0Muon_Cuts_Eta_5_10->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 9.0) && (tkRef->pt() <= 11.0) ) UnmatchedME0Muon_Cuts_Eta_9_11->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 10.0) && (tkRef->pt() <= 20.0) ) UnmatchedME0Muon_Cuts_Eta_10_50->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 20.0) && (tkRef->pt() <= 40.0) ) UnmatchedME0Muon_Cuts_Eta_50_100->Fill(fabs(tkRef->eta())); + if ( tkRef->pt() > 40.0) UnmatchedME0Muon_Cuts_Eta_100->Fill(fabs(tkRef->eta())); + + UnmatchedME0Muon_Cuts_Eta->Fill(fabs(tkRef->eta())); + UnmatchedME0Muon_Cuts_WideBinning_Eta->Fill(fabs(tkRef->eta())); + UnmatchedME0Muon_Cuts_WidestBinning_Eta->Fill(fabs(tkRef->eta())); + + UnmatchedME0Muon_ScatterPlot->Fill(fabs(tkRef->eta()), fabs(tkRef->phi()) ); + + for(unsigned int i=0; ieta()) < 2.8) ) UnmatchedME0Muon_Window_Pt->Fill(tkRef->pt()); + } + } + if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 2.8) ) { + UnmatchedME0Muon_Pt->Fill(tkRef->pt()); + UnmatchedME0Muon_SmallBins_Pt->Fill(tkRef->pt()); + UnmatchedME0Muon_VariableBins_Pt->Fill(tkRef->pt()); + } + + } + //if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 3.4) ) UnmatchedME0Muon_Pt->Fill(tkRef->pt()); + + } + ME0MuonID++; + } + + + + //Track Association by Chi2: + + //Map the list of all me0muons that failed or passed delR matching to a list of only Tight me0Muons that failed or passed delR matching + std::vector SkimmedIsMatched; + int i_me0muon=0; + for (std::vector::const_iterator thisMuon = OurMuons->begin(); + thisMuon != OurMuons->end(); ++thisMuon, ++i_me0muon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; + SkimmedIsMatched.push_back(IsMatched[i_me0muon]); + } //int w=0; //std::cout<<"associators size = "<(associator[ww]); + //associatorByChi2 = dynamic_cast(associator[ww]); + + + //associatorByHits = dynamic_cast(associator[ww]); //associatorByChi2 = associator[ww]; //std::cout<<"here now"<associateGenToReco(trackCollection, + // genParticles, + // &iEvent, + // &iSetup); + + // recSimColl=associatorByChi2->associateRecoToSim(trackCollection, + // trackingParticles, + // &iEvent, + // &iSetup); + + // simRecColl=associatorByChi2->associateSimToReco(trackCollection, + // trackingParticles, + // &iEvent, + // &iSetup); + + + recSimColl=associator[ww]->associateRecoToSim(trackCollection, + trackingParticles, &iEvent, &iSetup); - //std::cout<<"here now"<associateGenToReco(trackCollection, - genParticles, + + simRecColl=associator[ww]->associateSimToReco(trackCollection, + trackingParticles, &iEvent, &iSetup); + + // recSimColl=associatorByHits->associateRecoToSim(trackCollection, + // trackingParticles, + // &iEvent, + // &iSetup); + + // simRecColl=associatorByHits->associateSimToReco(trackCollection, + // trackingParticles, + // &iEvent, + // &iSetup); } //std::cout<<"here now"<status()<genParticle_begin())->numberOfMothers()<status() !=-99){ + int motherid=-1; + if ((*tp->genParticle_begin())->numberOfMothers()>0) { + if ((*tp->genParticle_begin())->mother()->numberOfMothers()>0){ + motherid=(*tp->genParticle_begin())->mother()->mother()->pdgId(); + } + } + + std::cout<<"Mother ID = "<status()==1) && ( (*tp->genParticle_begin())->numberOfMothers()==0 ) ) || + ( (tp->status()==1) ) ) SignalMuon=true; + + //if (SignalMuon) PreMatch_TP_R->Fill( sqrt(pow(tp->vertex().x(),2) + pow(tp->vertex().y(),2)) ); + } + //SignalMuon=true; + //if (tp->NumberOfMothers()==0) SignalMuon=true; + + if (SignalMuon){ + TPMuon_Eta->Fill(fabs(tp->eta())); + TPMuon_WideBinning_Eta->Fill(fabs(tp->eta())); + TPMuon_WidestBinning_Eta->Fill(fabs(tp->eta())); + if ( (fabs(tp->eta()) > 2.0) && (fabs(tp->eta()) < 2.8) ) { + TPMuon_SmallBins_Pt->Fill(tp->pt()); + TPMuon_VariableBins_Pt->Fill(tp->pt()); + TPMuon_Pt->Fill(tp->pt()); + } + if ( (tp->pt() > 5.0) && (tp->pt() <= 10.0) ) TPMuon_Eta_5_10->Fill(fabs(tpr->eta())); + if ( (tp->pt() > 10.0) && (tp->pt() <= 20.0) ) TPMuon_Eta_10_20->Fill(fabs(tpr->eta())); + if ( (tp->pt() > 20.0) && (tp->pt() <= 40.0) ) TPMuon_Eta_20_40->Fill(fabs(tpr->eta())); + if ( tp->pt() > 40.0) TPMuon_Eta_40->Fill(fabs(tpr->eta())); + + + } + + } + //std::cout<<"Section complete"<eta()) > 2.0) && (fabs(tp->eta()) < 2.8) ) PreMatch_TP_R->Fill( sqrt(pow(tp->vertex().x(),2) + pow(tp->vertex().y(),2)) ); + std::vector, double> > rt; //Check if the gen particle has been associated to any reco track @@ -771,12 +1537,12 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup //We take the first element of the vector, .begin(), and the trackRef from it, ->first, this is our best possible track match RefToBase assoc_recoTrack = rt.begin()->first; - std::cout<<"-----------------------------associated Track #"<second; - std::cout << "TrackingParticle #" <second; + //std::cout << "TrackingParticle #" < FakeRatePtCut) && (TMath::Abs(tp->eta()) < 2.8) )Chi2MatchedME0Muon_Eta->Fill(tp->eta()); @@ -785,58 +1551,168 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup } }//END if(simRecColl.find(tpr) != simRecColl.end()) - }//END for (GenParticleCollection::size_type i=0; isize(); i++){ + //END for (GenParticleCollection::size_type i=0; i::size_type i=0; i track(trackCollection, i); //std::vector > tp; - std::vector > tp; - std::vector > tpforfake; + //std::vector > tp; + //std::vector > tpforfake; + + std::vector > tp; + std::vector > tpforfake; //TrackingParticleRef tpr; - GenParticleRef tpr; - GenParticleRef tprforfake; + //GenParticleRef tpr; + //GenParticleRef tprforfake; + TrackingParticleRef tpr; + TrackingParticleRef tprforfake; //Check if the track is associated to any gen particle + bool TrackIsEfficient = false; + //std::cout<<"About to check first collection"<first; + //if (abs(tpr->pdgId()) != 13) continue; double assocChi2 = -(tp.begin()->second); //So this track is matched to a gen particle, lets get that gen particle now + if ( (simRecColl.find(tpr) != simRecColl.end()) ){ std::vector, double> > rt; - std::cout<<"Comparing gen and reco tracks"<Fill(fabs(tpr->eta())); + //std::cout<<"Comparing gen and reco tracks"< 0){ rt=simRecColl[tpr]; RefToBase bestrecotrackforeff = rt.begin()->first; //Only fill the efficiency histo if the track found matches up to a gen particle's best choice - if (bestrecotrackforeff == track) { - if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) )Chi2MatchedME0Muon_Eta->Fill(fabs(tpr->eta())); + if ( (bestrecotrackforeff == track ) && (abs(tpr->pdgId()) == 13) ) { + TrackIsEfficient=true; + //This section fills the numerator of the efficiency calculation... + //if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) + PostMatch_TP_R->Fill( sqrt(pow(tpr->vertex().x(),2) + pow(tpr->vertex().y(),2)) ); + if (tpr->eventId().bunchCrossing()) PostMatch_BX0_TP_R->Fill( sqrt(pow(tpr->vertex().x(),2) + pow(tpr->vertex().y(),2)) ); + + + if ( (tpr->pt() > FakeRatePtCut) ) + { + + + bool SignalMuon=false; + // for(unsigned int i=0; ipt()) < 0.25) && (CurrentParticle.numberOfMothers()==0) && (tpr->status()==1)) SignalMuon=true; + // } + + //if ( (tpr->status()==1) && ( tpr->genParticle_begin()->numberOfMothers()==0 ) ) SignalMuon=true; + //if ( (tpr->status()==1) && ( (*tpr->genParticle_begin())->numberOfMothers()==0 ) ) SignalMuon=true; + + //int motherid=(*tpr->genParticle_begin())->mother()->pdgId(); + + if (tpr->status() !=-99){ + int motherid=-1; + if ((*tpr->genParticle_begin())->numberOfMothers()>0) { + if ((*tpr->genParticle_begin())->mother()->numberOfMothers()>0){ + motherid=(*tpr->genParticle_begin())->mother()->mother()->pdgId(); + } + } + std::cout<<"Mother ID = "<status()==1) && ( (*tpr->genParticle_begin())->numberOfMothers()==0 ) ) || + ( (tpr->status()==1) ) )SignalMuon=true; + //( (tpr->status()==1) && ( (fabs(motherid)==22) || (fabs(motherid)==23) ) ) ) SignalMuon=true; + + } + //SignalMuon=true; + + if (SignalMuon){ + Chi2MatchedME0Muon_Eta->Fill(fabs(tpr->eta())); + Chi2MatchedME0Muon_WideBinning_Eta->Fill(fabs(tpr->eta())); + Chi2MatchedME0Muon_WidestBinning_Eta->Fill(fabs(tpr->eta())); + if ( (TMath::Abs(tpr->eta()) > 2.0) && (TMath::Abs(tpr->eta()) < 2.8) ) { + Chi2MatchedME0Muon_Pt->Fill(tpr->pt()); + Chi2MatchedME0Muon_SmallBins_Pt->Fill(tpr->pt()); + Chi2MatchedME0Muon_VariableBins_Pt->Fill(tpr->pt()); + } + + std::cout<<"eta = "<eta()<pt() > 5.0) && (track->pt() <= 10.0) ) Chi2MatchedME0Muon_Eta_5_10->Fill(fabs(tpr->eta())); + if ( (track->pt() > 9.0) && (track->pt() <= 11.0) ) Chi2MatchedME0Muon_Eta_9_11->Fill(fabs(tpr->eta())); + if ( (track->pt() > 10.0) && (track->pt() <= 20.0) ) Chi2MatchedME0Muon_Eta_10_20->Fill(fabs(tpr->eta())); + if ( (track->pt() > 20.0) && (track->pt() <= 40.0) ) Chi2MatchedME0Muon_Eta_20_40->Fill(fabs(tpr->eta())); + if ( track->pt() > 40.0) Chi2MatchedME0Muon_Eta_40->Fill(fabs(tpr->eta())); + } + + } + //...end section + if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) )AssociatedChi2_h->Fill(assocChi2); if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) )AssociatedChi2_Prob_h->Fill(TMath::Prob((assocChi2)*5,5)); - std::cout<<"assocChi2 = "<pt() >FakeRatePtCut) { + // if (tp.size()!=0) std::cout<<"Failed to find matching sim track: "<pt()<<", "<eta()<<", "<second<pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) { + Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + Chi2UnmatchedME0Muon_WideBinning_Eta->Fill(fabs(track->eta())); + Chi2UnmatchedME0Muon_WidestBinning_Eta->Fill(fabs(track->eta())); + if ( (TMath::Abs(track->eta()) > 2.0) && (TMath::Abs(track->eta()) < 2.8) ) { + Chi2UnmatchedME0Muon_Pt->Fill(track->pt()); + Chi2UnmatchedME0Muon_SmallBins_Pt->Fill(track->pt()); + Chi2UnmatchedME0Muon_VariableBins_Pt->Fill(track->pt()); + } } + //std::cout<<"unmatched eta = "<eta(); + } //End checking of Efficient muons - //For Fakes -------------------------------------------- + //For Fakes -------------------------------------------- here we fill the numerator for the F.R., Chi2UnmatchedME0Muon_Eta + //The denominator is filled elsewhere, just a histo of all the ME0Muon eta values + //It is ME0Muon_Cuts_Eta, so named because it is all ME0Muons passing the selection (also the pT cut) //Check if finding a track associated to a gen particle fails, or if there is no track in the collection at all if( (recSimColl.find(track) == recSimColl.end() ) || ( recSimColl[track].size() == 0 ) ){ //So we've now determined that this track is not associated to any gen, and fill our histo of fakes: - if ((track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // if ((track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) { + // Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // std::cout<<"unmatched eta = "<eta(); + // } + //Check if this muon matched via Del-R matching + if (SkimmedIsMatched[i]){ + + + if ((track->pt() >FakeRatePtCut) ){ + if (tp.size()!=0) std::cout<<"Found an me0muontrack failing chi2 matching: "<pt()<<", "<eta()<<", "<second<first; //We now have the gen particle, to check + //If for some crazy reason we can't find the gen particle, that means its a fake if ( (simRecColl.find(tprforfake) == simRecColl.end()) || (simRecColl[tprforfake].size() == 0) ) { - if ((track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8)) Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // if ((track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8)) { + // Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // std::cout<<"unmatched eta = "<eta(); + // } + //Check if this muon matched via Del-R matching + if (SkimmedIsMatched[i]){ + if ((track->pt() >FakeRatePtCut) ) { + if (tp.size()!=0) std::cout<<"Found an me0muontrack failing chi2 matching: "<pt()<<", "<eta()<<", "<second< 0) { @@ -859,7 +1745,14 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup RefToBase bestrecotrack = rtforfake.begin()->first; //if the best reco track is NOT the track that we're looking at, we know we have a fake, that was within the cut, but not the closest if (bestrecotrack != track) { - if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // if ( (track->pt() > FakeRatePtCut) && (TMath::Abs(track->eta()) < 2.8) ) { + // Chi2UnmatchedME0Muon_Eta->Fill(fabs(track->eta())); + // std::cout<<"unmatched eta = "<eta(); + // } + //Check if this muon matched via Del-R matching + if (IsMatched[i]){ + if (tp.size()!=0) std::cout<<"Found an me0muontrack failing chi2 matching: "<pt()<<", "<eta()<<", "<second<::size_type i=0; i MatchedSegIds; - - for(unsigned int i=0; i 2.0) && (fabs(CurrentParticle.eta()) < 2.8) ) GenMuon_Eta->Fill(CurrentParticle.eta()); - //std::cout<<"Mother's ID is: "< ReferenceTrackPt; - - double VertexDiff=-1,PtDiff=-1,QOverPtDiff=-1,PDiff=-1; - - std::cout<<"Size = "<size()<::const_iterator thisMuon = OurMuons->begin(); - thisMuon != OurMuons->end(); ++thisMuon){ - TrackRef tkRef = thisMuon->innerTrack(); - SegIdForMatch.push_back(thisMuon->me0segid()); - thisDelR = reco::deltaR(CurrentParticle,*tkRef); - ReferenceTrackPt.push_back(tkRef->pt()); - if (tkRef->pt() > FakeRatePtCut ) { - if (thisDelR < MatchingWindowDelR ){ - if (thisDelR < LowestDelR){ - LowestDelR = thisDelR; - //if (fabs(tkRef->pt() - CurrentParticle.pt())/CurrentParticle.pt() < 0.50) MatchedID = ME0MuonID; - MatchedID = ME0MuonID; - VertexDiff = fabs(tkRef->vz()-CurrentParticle.vz()); - QOverPtDiff = ( (tkRef->charge() /tkRef->pt()) - (CurrentParticle.charge()/CurrentParticle.pt() ) )/ (CurrentParticle.charge()/CurrentParticle.pt() ); - PtDiff = (tkRef->pt() - CurrentParticle.pt())/CurrentParticle.pt(); - PDiff = (tkRef->p() - CurrentParticle.p())/CurrentParticle.p(); - } - } - } - - ME0MuonID++; - - } - - if (MatchedID != -1){ - IsMatched[MatchedID] = true; - if (CurrentParticle.pt() >FakeRatePtCut) { - MatchedME0Muon_Eta->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) MatchedME0Muon_Eta_5_10->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) MatchedME0Muon_Eta_10_20->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) MatchedME0Muon_Eta_20_40->Fill(fabs(CurrentParticle.eta())); - if ( CurrentParticle.pt() > 40.0) MatchedME0Muon_Eta_40->Fill(fabs(CurrentParticle.eta())); - - - - VertexDiff_h->Fill(VertexDiff); - PtDiff_h->Fill(PtDiff); - QOverPtDiff_h->Fill(QOverPtDiff); - PtDiff_s->Fill(CurrentParticle.eta(),PtDiff); - if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) PtDiff_s_5_10->Fill(CurrentParticle.eta(),PtDiff); - if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) PtDiff_s_10_20->Fill(CurrentParticle.eta(),PtDiff); - if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) PtDiff_s_20_40->Fill(CurrentParticle.eta(),PtDiff); - if ( CurrentParticle.pt() > 40.0) PtDiff_s_40->Fill(CurrentParticle.eta(),PtDiff); - PtDiff_p->Fill(CurrentParticle.eta(),PtDiff); - - PDiff_h->Fill(PDiff); - PDiff_s->Fill(CurrentParticle.eta(),PDiff); - PDiff_p->Fill(CurrentParticle.eta(),PDiff); - } - MatchedSegIds.push_back(SegIdForMatch[MatchedID]); - - // if ( ((CurrentParticle.eta()) > 2.4) && ((CurrentParticle.eta()) < 3.8) ) { - // //MatchedME0Muon_Pt->Fill(CurrentParticle.pt()); - // //MatchedME0Muon_Pt->Fill(.pt()); - - // } - } - if (CurrentParticle.pt() >FakeRatePtCut) { - GenMuon_Eta->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 5.0) && (CurrentParticle.pt() <= 10.0) ) GenMuon_Eta_5_10->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 10.0) && (CurrentParticle.pt() <= 20.0) ) GenMuon_Eta_10_20->Fill(fabs(CurrentParticle.eta())); - if ( (CurrentParticle.pt() > 20.0) && (CurrentParticle.pt() <= 40.0) ) GenMuon_Eta_20_40->Fill(fabs(CurrentParticle.eta())); - if ( CurrentParticle.pt() > 40.0) GenMuon_Eta_40->Fill(fabs(CurrentParticle.eta())); - GenMuon_Phi->Fill(CurrentParticle.phi()); - if ( ((CurrentParticle.eta()) > 2.0) && ((CurrentParticle.eta()) < 2.8) ) GenMuon_Pt->Fill(CurrentParticle.pt()); - } - - } - } - - -//Diff studies for gen matching - - std::cout<<"Doing first propagation"< 2.8 )) continue; - - float zSign = CurrentParticle.pz()/fabs(CurrentParticle.pz()); - - float zValue = 526.75 * zSign; - Plane *plane = new Plane(Surface::PositionType(0,0,zValue),Surface::RotationType()); - TLorentzVector Momentum; - Momentum.SetPtEtaPhiM(CurrentParticle.pt() - ,CurrentParticle.eta() - ,CurrentParticle.phi() - ,CurrentParticle.mass()); - GlobalVector p3gen(Momentum.Px(), Momentum.Py(), Momentum.Pz()); - GlobalVector r3gen = GlobalVector(CurrentParticle.vertex().x() - ,CurrentParticle.vertex().y() - ,CurrentParticle.vertex().z()); - - AlgebraicSymMatrix66 covGen = AlgebraicMatrixID(); - covGen *= 1e-20; // initialize to sigma=1e-10 .. should get overwhelmed by MULS - AlgebraicSymMatrix66 covFinal; - int chargeGen = CurrentParticle.charge(); - - //Propagation - FreeTrajectoryState initstate = getFTS(p3gen, r3gen, chargeGen, covGen, &*bField); - - SteppingHelixStateInfo startstate(initstate); - SteppingHelixStateInfo laststate; - - const SteppingHelixPropagator* ThisshProp = - dynamic_cast(&*shProp); - - laststate = ThisshProp->propagate(startstate, *plane); - - FreeTrajectoryState finalstate; - laststate.getFreeState(finalstate); - - GlobalVector p3Final, r3Final; - getFromFTS(finalstate, p3Final, r3Final, chargeGen, covFinal); - - - int ME0MuonID = 0; - - for (std::vector::const_iterator thisMuon = OurMuons->begin(); - thisMuon != OurMuons->end(); ++thisMuon){ - - TrackRef tkRef = thisMuon->innerTrack(); - SegIdForMatch.push_back(thisMuon->me0segid()); - - ME0Segment Seg = thisMuon->me0segment(); - ME0DetId id =Seg.me0DetId(); - auto roll = me0Geom->etaPartition(id); - - double DirectionPull, DirectionPullNum, DirectionPullDenom; - - //Computing the sigma for the track direction - Double_t mag_track = p3Final.perp(); - //Double_t phi_track = p3Final.phi(); - - //Double_t dmagdx_track = p3Final.x()/mag_track; - //Double_t dmagdy_track = p3Final.y()/mag_track; - Double_t dphidx_track = -p3Final.y()/(mag_track*mag_track); - Double_t dphidy_track = p3Final.x()/(mag_track*mag_track); - Double_t sigmaphi_track = sqrt( dphidx_track*dphidx_track*covFinal(3,3)+ - dphidy_track*dphidy_track*covFinal(4,4)+ - dphidx_track*dphidy_track*2*covFinal(3,4) ); - - DirectionPullNum = p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi(); - DirectionPullDenom = sqrt( pow(roll->toGlobal(Seg.localPosition()).phi(),2) + pow(sigmaphi_track,2) ); - DirectionPull = DirectionPullNum / DirectionPullDenom; - - - if (IsMatched[ME0MuonID]){ - SegGenDirPhiDiff_True_h->Fill(p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegGenDirEtaDiff_True_h->Fill(p3Final.eta()-roll->toGlobal(Seg.localDirection()).eta() ); - SegGenDirPhiPull_True_h->Fill(DirectionPull); - } - - if ((zSign * roll->toGlobal(Seg.localDirection()).z()) > 0 ){ - SegGenDirPhiDiff_All_h->Fill(p3Final.phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegGenDirPhiPull_All_h->Fill(DirectionPull); - - SegGenDirEtaDiff_All_h->Fill(p3Final.eta()-roll->toGlobal(Seg.localDirection()).eta() ); - } - ME0MuonID++; - } - } - } - - //Del R study =========================== - for(unsigned int i=0; i::const_iterator thisMuon = OurMuons->begin(); - thisMuon != OurMuons->end(); ++thisMuon){ - TrackRef tkRef = thisMuon->innerTrack(); - thisDelR = reco::deltaR(CurrentParticle,*tkRef); - if (thisDelR < LowestDelR) LowestDelR = thisDelR; - } - - ClosestDelR_s->Fill(CurrentParticle.eta(), LowestDelR); - ClosestDelR_p->Fill(CurrentParticle.eta(), LowestDelR); - } - } - - //==================================== - - // -----Finally, we loop over all the ME0Muons in the event - // -----Before, we plotted the gen muon pt and eta for the efficiency plot of matches - // -----Now, each time a match failed, we plot the ME0Muon pt and eta - int ME0MuonID = 0; - for (std::vector::const_iterator thisMuon = OurMuons->begin(); - thisMuon != OurMuons->end(); ++thisMuon){ - TrackRef tkRef = thisMuon->innerTrack(); - - - if (IsMatched[ME0MuonID]) { - if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 2.8) ) MatchedME0Muon_Pt->Fill(tkRef->pt()); - - //Moved resolution stuff here, only calculate resolutions for matched muons! - - - } - - if (!IsMatched[ME0MuonID]){ - - UnmatchedME0Muon_Eta->Fill(fabs(tkRef->eta())); - if ((tkRef->pt() > FakeRatePtCut) && (TMath::Abs(tkRef->eta()) < 2.8) ) { - if ( (tkRef->pt() > 5.0) && (tkRef->pt() <= 10.0) ) UnmatchedME0Muon_Cuts_Eta_5_10->Fill(fabs(tkRef->eta())); - if ( (tkRef->pt() > 10.0) && (tkRef->pt() <= 20.0) ) UnmatchedME0Muon_Cuts_Eta_10_20->Fill(fabs(tkRef->eta())); - if ( (tkRef->pt() > 20.0) && (tkRef->pt() <= 40.0) ) UnmatchedME0Muon_Cuts_Eta_20_40->Fill(fabs(tkRef->eta())); - if ( tkRef->pt() > 40.0) UnmatchedME0Muon_Cuts_Eta_40->Fill(fabs(tkRef->eta())); - - UnmatchedME0Muon_Cuts_Eta->Fill(fabs(tkRef->eta())); - } - //if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 3.4) ) UnmatchedME0Muon_Pt->Fill(tkRef->pt()); - if ( (TMath::Abs(tkRef->eta()) < 2.8) ) UnmatchedME0Muon_Pt->Fill(tkRef->pt()); - } - ME0MuonID++; - } - - + //delete recSimColl; + //delete simRecColl; + }//END for (unsigned int www=0;www::const_iterator thisSegment = OurSegments->begin(); // thisSegment != OurSegments->end();++thisSegment){ @@ -1142,6 +1798,59 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup } + // //Standalone Muon study: + + // std::cout<<"About to get muons:"< > muons; + // iEvent.getByLabel("Muons", muons); + + // std::cout<<"Have muons, about to start"< ReferenceTrackPt; + + // //double VertexDiff=-1,PtDiff=-1,QOverPtDiff=-1,PDiff=-1,MatchedEta=-1; + // double PtDiff=-1,QOverPtDiff=-1,MatchedEta=-1; + + // std::cout<<"Size = "<size()<::const_iterator thisMuon = muons->begin(); + // thisMuon != muons->end(); ++thisMuon){ + // if (thisMuon->isStandAloneMuon()){ + + // TrackRef tkRef = thisMuon->outerTrack(); + // thisDelR = reco::deltaR(CurrentParticle,*tkRef); + // if (tkRef->pt() > FakeRatePtCut ) { + // if (thisDelR < MatchingWindowDelR ){ + // if (thisDelR < LowestDelR){ + // LowestDelR = thisDelR; + + // MatchedEta=fabs(CurrentParticle.eta()); + // //VertexDiff = fabs(tkRef->vz()-CurrentParticle.vz()); + // QOverPtDiff = ( (tkRef->charge() /tkRef->pt()) - (CurrentParticle.charge()/CurrentParticle.pt() ) )/ (CurrentParticle.charge()/CurrentParticle.pt() ); + // PtDiff = (tkRef->pt() - CurrentParticle.pt())/CurrentParticle.pt(); + // //PDiff = (tkRef->p() - CurrentParticle.p())/CurrentParticle.p(); + // } + // } + // } + // } + // } + // StandaloneMatchedME0Muon_Eta->Fill(MatchedEta); + // StandaloneMatchedME0Muon_WideBinning_Eta->Fill(MatchedEta); + // StandaloneMatchedME0Muon_WidestBinning_Eta->Fill(MatchedEta); + + // //StandaloneVertexDiff_h->Fill(VertexDiff); + // StandalonePtDiff_h->Fill(PtDiff); + // StandaloneQOverPtDiff_h->Fill(QOverPtDiff); + // StandalonePtDiff_s->Fill(CurrentParticle.eta(),PtDiff); + + // } + // } + // for (std::vector::const_iterator thisMuon = OurMuons->begin(); // thisMuon != OurMuons->end(); ++thisMuon){ // TrackRef tkRef = thisMuon->innerTrack(); @@ -1158,10 +1867,11 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup std::vector UniqueIdList; int TrackID=0; - std::cout<<"Doing some propagation"<::const_iterator thisMuon = OurMuons->begin(); thisMuon != OurMuons->end(); ++thisMuon){ + if (!muon::isGoodMuon(me0Geom, *thisMuon, muon::Tight)) continue; TrackRef tkRef = thisMuon->innerTrack(); //ME0Segment segRef = thisMuon->me0segment(); //const ME0Segment* SegId = segRef->get(); @@ -1177,114 +1887,114 @@ ME0MuonAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup //For a direction study... - if ( (tkRef->pt() > 5.0) && (IsMatched[MuID]) ){ - SegPosDirPhiDiff_True_h->Fill(roll->toGlobal(Seg.localPosition()).phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegPosDirEtaDiff_True_h->Fill(roll->toGlobal(Seg.localPosition()).eta()-roll->toGlobal(Seg.localDirection()).eta() ); - } + // if ( (tkRef->pt() > 5.0) && (IsMatched[MuID]) ){ + // SegPosDirPhiDiff_True_h->Fill(roll->toGlobal(Seg.localPosition()).phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegPosDirEtaDiff_True_h->Fill(roll->toGlobal(Seg.localPosition()).eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // } - SegPosDirPhiDiff_All_h->Fill(roll->toGlobal(Seg.localPosition()).phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegPosDirEtaDiff_All_h->Fill(roll->toGlobal(Seg.localPosition()).eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // SegPosDirPhiDiff_All_h->Fill(roll->toGlobal(Seg.localPosition()).phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegPosDirEtaDiff_All_h->Fill(roll->toGlobal(Seg.localPosition()).eta()-roll->toGlobal(Seg.localDirection()).eta() ); - // For another direction study... - float zSign = tkRef->pz()/fabs(tkRef->pz()); + // // For another direction study... + // float zSign = tkRef->pz()/fabs(tkRef->pz()); - //float zValue = 560. * zSign; - float zValue = 526.75 * zSign; - Plane *plane = new Plane(Surface::PositionType(0,0,zValue),Surface::RotationType()); - //Getting the initial variables for propagation - int chargeReco = tkRef->charge(); - GlobalVector p3reco, r3reco; + // //float zValue = 560. * zSign; + // float zValue = 526.75 * zSign; + // Plane *plane = new Plane(Surface::PositionType(0,0,zValue),Surface::RotationType()); + // //Getting the initial variables for propagation + // int chargeReco = tkRef->charge(); + // GlobalVector p3reco, r3reco; - p3reco = GlobalVector(tkRef->outerPx(), tkRef->outerPy(), tkRef->outerPz()); - r3reco = GlobalVector(tkRef->outerX(), tkRef->outerY(), tkRef->outerZ()); + // p3reco = GlobalVector(tkRef->outerPx(), tkRef->outerPy(), tkRef->outerPz()); + // r3reco = GlobalVector(tkRef->outerX(), tkRef->outerY(), tkRef->outerZ()); - AlgebraicSymMatrix66 covReco; - //This is to fill the cov matrix correctly - AlgebraicSymMatrix55 covReco_curv; - covReco_curv = tkRef->outerStateCovariance(); - FreeTrajectoryState initrecostate = getFTS(p3reco, r3reco, chargeReco, covReco_curv, &*bField); - getFromFTS(initrecostate, p3reco, r3reco, chargeReco, covReco); + // AlgebraicSymMatrix66 covReco; + // //This is to fill the cov matrix correctly + // AlgebraicSymMatrix55 covReco_curv; + // covReco_curv = tkRef->outerStateCovariance(); + // FreeTrajectoryState initrecostate = getFTS(p3reco, r3reco, chargeReco, covReco_curv, &*bField); + // getFromFTS(initrecostate, p3reco, r3reco, chargeReco, covReco); - //Now we propagate and get the propagated variables from the propagated state - SteppingHelixStateInfo startrecostate(initrecostate); - SteppingHelixStateInfo lastrecostate; + // //Now we propagate and get the propagated variables from the propagated state + // SteppingHelixStateInfo startrecostate(initrecostate); + // SteppingHelixStateInfo lastrecostate; - const SteppingHelixPropagator* ThisshProp = - dynamic_cast(&*shProp); + // const SteppingHelixPropagator* ThisshProp = + // dynamic_cast(&*shProp); - lastrecostate = ThisshProp->propagate(startrecostate, *plane); + // lastrecostate = ThisshProp->propagate(startrecostate, *plane); - FreeTrajectoryState finalrecostate; - lastrecostate.getFreeState(finalrecostate); + // FreeTrajectoryState finalrecostate; + // lastrecostate.getFreeState(finalrecostate); - AlgebraicSymMatrix66 covFinalReco; - GlobalVector p3FinalReco_glob, r3FinalReco_globv; - getFromFTS(finalrecostate, p3FinalReco_glob, r3FinalReco_globv, chargeReco, covFinalReco); - GlobalPoint r3FinalReco_glob(r3FinalReco_globv.x(),r3FinalReco_globv.y(),r3FinalReco_globv.z()); - - double DirectionPull, DirectionPullNum, DirectionPullDenom; - - //Computing the sigma for the track direction - Double_t mag_track = p3FinalReco_glob.perp(); - //Double_t phi_track = p3FinalReco_glob.phi(); - - //Double_t dmagdx_track = p3FinalReco_glob.x()/mag_track; - //Double_t dmagdy_track = p3FinalReco_glob.y()/mag_track; - Double_t dphidx_track = -p3FinalReco_glob.y()/(mag_track*mag_track); - Double_t dphidy_track = p3FinalReco_glob.x()/(mag_track*mag_track); - Double_t sigmaphi_track = sqrt( dphidx_track*dphidx_track*covFinalReco(3,3)+ - dphidy_track*dphidy_track*covFinalReco(4,4)+ - dphidx_track*dphidy_track*2*covFinalReco(3,4) ); - - DirectionPullNum = p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi(); - DirectionPullDenom = sqrt( pow(roll->toGlobal(Seg.localPosition()).phi(),2) + pow(sigmaphi_track,2) ); - DirectionPull = DirectionPullNum / DirectionPullDenom; - - if ( (tkRef->pt() > 5.0)&& (IsMatched[MuID]) ){ - SegTrackDirPhiDiff_True_h->Fill(p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegTrackDirEtaDiff_True_h->Fill(p3FinalReco_glob.eta()-roll->toGlobal(Seg.localDirection()).eta() ); - SegTrackDirPhiPull_True_h->Fill(DirectionPull); - } - SegTrackDirPhiDiff_All_h->Fill(p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi() ); - SegTrackDirPhiPull_All_h->Fill(DirectionPull); + // AlgebraicSymMatrix66 covFinalReco; + // GlobalVector p3FinalReco_glob, r3FinalReco_globv; + // getFromFTS(finalrecostate, p3FinalReco_glob, r3FinalReco_globv, chargeReco, covFinalReco); + // GlobalPoint r3FinalReco_glob(r3FinalReco_globv.x(),r3FinalReco_globv.y(),r3FinalReco_globv.z()); + + // double DirectionPull, DirectionPullNum, DirectionPullDenom; + + // //Computing the sigma for the track direction + // Double_t mag_track = p3FinalReco_glob.perp(); + // //Double_t phi_track = p3FinalReco_glob.phi(); + + // //Double_t dmagdx_track = p3FinalReco_glob.x()/mag_track; + // //Double_t dmagdy_track = p3FinalReco_glob.y()/mag_track; + // Double_t dphidx_track = -p3FinalReco_glob.y()/(mag_track*mag_track); + // Double_t dphidy_track = p3FinalReco_glob.x()/(mag_track*mag_track); + // Double_t sigmaphi_track = sqrt( dphidx_track*dphidx_track*covFinalReco(3,3)+ + // dphidy_track*dphidy_track*covFinalReco(4,4)+ + // dphidx_track*dphidy_track*2*covFinalReco(3,4) ); + + // DirectionPullNum = p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi(); + // DirectionPullDenom = sqrt( pow(roll->toGlobal(Seg.localPosition()).phi(),2) + pow(sigmaphi_track,2) ); + // DirectionPull = DirectionPullNum / DirectionPullDenom; + + // if ( (tkRef->pt() > 5.0)&& (IsMatched[MuID]) ){ + // SegTrackDirPhiDiff_True_h->Fill(p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegTrackDirEtaDiff_True_h->Fill(p3FinalReco_glob.eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // SegTrackDirPhiPull_True_h->Fill(DirectionPull); + // } + // SegTrackDirPhiDiff_All_h->Fill(p3FinalReco_glob.phi()-roll->toGlobal(Seg.localDirection()).phi() ); + // SegTrackDirPhiPull_All_h->Fill(DirectionPull); - SegTrackDirEtaDiff_All_h->Fill(p3FinalReco_glob.eta()-roll->toGlobal(Seg.localDirection()).eta() ); + // SegTrackDirEtaDiff_All_h->Fill(p3FinalReco_glob.eta()-roll->toGlobal(Seg.localDirection()).eta() ); - LocalPoint r3FinalReco = roll->toLocal(r3FinalReco_glob); - LocalVector p3FinalReco=roll->toLocal(p3FinalReco_glob); - LocalTrajectoryParameters ltp(r3FinalReco,p3FinalReco,chargeReco); - JacobianCartesianToLocal jctl(roll->surface(),ltp); - AlgebraicMatrix56 jacobGlbToLoc = jctl.jacobian(); + // LocalPoint r3FinalReco = roll->toLocal(r3FinalReco_glob); + // LocalVector p3FinalReco=roll->toLocal(p3FinalReco_glob); + // LocalTrajectoryParameters ltp(r3FinalReco,p3FinalReco,chargeReco); + // JacobianCartesianToLocal jctl(roll->surface(),ltp); + // AlgebraicMatrix56 jacobGlbToLoc = jctl.jacobian(); - AlgebraicMatrix55 Ctmp = (jacobGlbToLoc * covFinalReco) * ROOT::Math::Transpose(jacobGlbToLoc); - AlgebraicSymMatrix55 C; // I couldn't find any other way, so I resort to the brute force - for(int i=0; i<5; ++i) { - for(int j=0; j<5; ++j) { - C[i][j] = Ctmp[i][j]; + // AlgebraicMatrix55 Ctmp = (jacobGlbToLoc * covFinalReco) * ROOT::Math::Transpose(jacobGlbToLoc); + // AlgebraicSymMatrix55 C; // I couldn't find any other way, so I resort to the brute force + // for(int i=0; i<5; ++i) { + // for(int j=0; j<5; ++j) { + // C[i][j] = Ctmp[i][j]; - } - } + // } + // } - LocalPoint thisPosition(Seg.localPosition()); + // LocalPoint thisPosition(Seg.localPosition()); - Double_t sigmax = sqrt(C[3][3]+Seg.localPositionError().xx() ); - Double_t sigmay = sqrt(C[4][4]+Seg.localPositionError().yy() ); + // Double_t sigmax = sqrt(C[3][3]+Seg.localPositionError().xx() ); + // Double_t sigmay = sqrt(C[4][4]+Seg.localPositionError().yy() ); - XPull_h->Fill((thisPosition.x()-r3FinalReco.x())/sigmax); - YPull_h->Fill((thisPosition.y()-r3FinalReco.y())/sigmay); + // XPull_h->Fill((thisPosition.x()-r3FinalReco.x())/sigmax); + // YPull_h->Fill((thisPosition.y()-r3FinalReco.y())/sigmay); - XDiff_h->Fill((thisPosition.x()-r3FinalReco.x())); - YDiff_h->Fill((thisPosition.y()-r3FinalReco.y())); + // XDiff_h->Fill((thisPosition.x()-r3FinalReco.x())); + // YDiff_h->Fill((thisPosition.y()-r3FinalReco.y())); - //std::cout<<"AM HERE"<pt() > FakeRatePtCut)&& (IsMatched[MuID]) ){ + // //std::cout<<"AM HERE"<pt() > FakeRatePtCut)&& (IsMatched[MuID]) ){ - //std::cout<<"thisPosition = "<pt() > FakeRatePtCut) && (TMath::Abs(tkRef->eta()) < 2.8)){ ME0Muon_Cuts_Eta->Fill(fabs(tkRef->eta())); + ME0Muon_Cuts_WideBinning_Eta->Fill(fabs(tkRef->eta())); + ME0Muon_Cuts_WidestBinning_Eta->Fill(fabs(tkRef->eta())); if ( (tkRef->pt() > 5.0) && (tkRef->pt() <= 10.0) ) ME0Muon_Cuts_Eta_5_10->Fill(fabs(tkRef->eta())); - if ( (tkRef->pt() > 10.0) && (tkRef->pt() <= 20.0) ) ME0Muon_Cuts_Eta_10_20->Fill(fabs(tkRef->eta())); - if ( (tkRef->pt() > 20.0) && (tkRef->pt() <= 40.0) ) ME0Muon_Cuts_Eta_20_40->Fill(fabs(tkRef->eta())); - if ( tkRef->pt() > 40.0) ME0Muon_Cuts_Eta_40->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 9.0) && (tkRef->pt() <= 11.0) ) ME0Muon_Cuts_Eta_9_11->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 10.0) && (tkRef->pt() <= 20.0) ) ME0Muon_Cuts_Eta_10_50->Fill(fabs(tkRef->eta())); + if ( (tkRef->pt() > 20.0) && (tkRef->pt() <= 40.0) ) ME0Muon_Cuts_Eta_50_100->Fill(fabs(tkRef->eta())); + if ( tkRef->pt() > 40.0) ME0Muon_Cuts_Eta_100->Fill(fabs(tkRef->eta())); + + + if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 2.8) ) { + ME0Muon_Pt->Fill(tkRef->pt()); + ME0Muon_SmallBins_Pt->Fill(tkRef->pt()); + ME0Muon_VariableBins_Pt->Fill(tkRef->pt()); + } } - if ( (TMath::Abs(tkRef->eta()) > 2.0) && (TMath::Abs(tkRef->eta()) < 2.8) ) ME0Muon_Pt->Fill(tkRef->pt()); - TrackID++; MuID++; - } + + } - std::cout<cd(); + TCanvas *c1 = new TCanvas("c1", "canvas" ); @@ -1503,6 +2226,7 @@ void ME0MuonAnalyzer::endJob() c1->Print("DiffY.png"); gStyle->SetOptStat(0); + Nevents_h->Write(); SegPosDirPhiDiff_True_h->Write(); SegPosDirPhiDiff_True_h->Draw(); c1->Print(histoFolder+"/SegPosDirPhiDiff_True_h.png"); SegPosDirEtaDiff_True_h->Write(); SegPosDirEtaDiff_True_h->Draw(); c1->Print(histoFolder+"/SegPosDirEtaDiff_True_h.png"); @@ -1565,37 +2289,108 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/CheckME0Muon_Eta.png"); ME0Muon_Cuts_Eta->Write(); ME0Muon_Cuts_Eta->Draw(); c1->Print(histoFolder+"/ME0Muon_Cuts_Eta.png"); + ME0Muon_Cuts_WidestBinning_Eta->Write(); ME0Muon_Cuts_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/ME0Muon_Cuts_WidestBinning_Eta.png"); + ME0Muon_Cuts_WideBinning_Eta->Write(); ME0Muon_Cuts_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/ME0Muon_Cuts_WideBinning_Eta.png"); //c1->SetLogy(); ME0Muon_Pt->Write(); ME0Muon_Pt->Draw(); ME0Muon_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); ME0Muon_Pt->GetXaxis()->SetTitleSize(0.05); c1->Print(histoFolder+"/ME0Muon_Pt.png"); + ME0Muon_SmallBins_Pt->Write(); ME0Muon_SmallBins_Pt->Draw(); + ME0Muon_SmallBins_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); + ME0Muon_SmallBins_Pt->GetXaxis()->SetTitleSize(0.05); + c1->Print(histoFolder+"/ME0Muon_SmallBins_Pt.png"); + + ME0Muon_VariableBins_Pt->Write(); ME0Muon_VariableBins_Pt->Draw(); + ME0Muon_VariableBins_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); + ME0Muon_VariableBins_Pt->GetXaxis()->SetTitleSize(0.05); + c1->Print(histoFolder+"/ME0Muon_VariableBins_Pt.png"); + GenMuon_Eta->Write(); GenMuon_Eta->Draw(); c1->Print(histoFolder+"/GenMuon_Eta.png"); + GenMuon_WideBinning_Eta->Write(); GenMuon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/GenMuon_WideBinning_Eta.png"); + GenMuon_WidestBinning_Eta->Write(); GenMuon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/GenMuon_WidestBinning_Eta.png"); + + TPMuon_Eta->Write(); TPMuon_Eta->Draw(); c1->Print(histoFolder+"/TPMuon_Eta.png"); + TPMuon_WideBinning_Eta->Write(); TPMuon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/TPMuon_WideBinning_Eta.png"); + TPMuon_WidestBinning_Eta->Write(); TPMuon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/TPMuon_WidestBinning_Eta.png"); + + + TPMuon_Pt->Write(); TPMuon_Pt->Draw(); c1->Print(histoFolder+"/TPMuon_Pt.png"); + TPMuon_SmallBins_Pt->Write(); TPMuon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/TPMuon_SmallBins_Pt.png"); + TPMuon_VariableBins_Pt->Write(); TPMuon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/TPMuon_VariableBins_Pt.png"); GenMuon_Pt->Write(); GenMuon_Pt->Draw(); c1->Print(histoFolder+"/GenMuon_Pt.png"); + GenMuon_SmallBins_Pt->Write(); GenMuon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/GenMuon_SmallBins_Pt.png"); + GenMuon_VariableBins_Pt->Write(); GenMuon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/GenMuon_VariableBins_Pt.png"); MatchedME0Muon_Eta->Write(); MatchedME0Muon_Eta->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_Eta.png"); + MatchedME0Muon_WideBinning_Eta->Write(); MatchedME0Muon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_WideBinning_Eta.png"); + MatchedME0Muon_WidestBinning_Eta->Write(); MatchedME0Muon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_WidestBinning_Eta.png"); + + StandaloneMatchedME0Muon_Eta->Write(); StandaloneMatchedME0Muon_Eta->Draw(); c1->Print(histoFolder+"/StandaloneMatchedME0Muon_Eta.png"); + StandaloneMatchedME0Muon_WideBinning_Eta->Write(); StandaloneMatchedME0Muon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/StandaloneMatchedME0Muon_WideBinning_Eta.png"); + StandaloneMatchedME0Muon_WidestBinning_Eta->Write(); StandaloneMatchedME0Muon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/StandaloneMatchedME0Muon_WidestBinning_Eta.png"); Chi2MatchedME0Muon_Eta->Write(); Chi2MatchedME0Muon_Eta->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_Eta.png"); + Chi2MatchedME0Muon_WideBinning_Eta->Write(); Chi2MatchedME0Muon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_WideBinning_Eta.png"); + Chi2MatchedME0Muon_WidestBinning_Eta->Write(); Chi2MatchedME0Muon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_WidestBinning_Eta.png"); Chi2UnmatchedME0Muon_Eta->Write(); Chi2UnmatchedME0Muon_Eta->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_Eta.png"); + Chi2UnmatchedME0Muon_WideBinning_Eta->Write(); Chi2UnmatchedME0Muon_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_WideBinning_Eta.png"); + Chi2UnmatchedME0Muon_WidestBinning_Eta->Write(); Chi2UnmatchedME0Muon_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_WidestBinning_Eta.png"); gStyle->SetOptStat(1); MatchedME0Muon_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); //MatchedME0Muon_Pt->GetYaxis()->SetTitle(" \# of Se"); MatchedME0Muon_Pt->Write(); MatchedME0Muon_Pt->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_Pt.png"); + + MatchedME0Muon_SmallBins_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); + //MatchedME0Muon_SmallBins_Pt->GetYaxis()->SetTitle(" \# of Se"); + + MatchedME0Muon_SmallBins_Pt->Write(); MatchedME0Muon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_SmallBins_Pt.png"); + + MatchedME0Muon_VariableBins_Pt->GetXaxis()->SetTitle("ME0Muon p_{T}"); + //MatchedME0Muon_VariableBins_Pt->GetYaxis()->SetTitle(" \# of Se"); + + MatchedME0Muon_VariableBins_Pt->Write(); MatchedME0Muon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/MatchedME0Muon_VariableBins_Pt.png"); gStyle->SetOptStat(0); UnmatchedME0Muon_Eta->Write(); UnmatchedME0Muon_Eta->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Eta.png"); UnmatchedME0Muon_Cuts_Eta->Write(); UnmatchedME0Muon_Cuts_Eta->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Cuts_Eta.png"); + UnmatchedME0Muon_Cuts_WideBinning_Eta->Write(); UnmatchedME0Muon_Cuts_WideBinning_Eta->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Cuts_WideBinning_Eta.png"); + UnmatchedME0Muon_Cuts_WidestBinning_Eta->Write(); UnmatchedME0Muon_Cuts_WidestBinning_Eta->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Cuts_WidestBinning_Eta.png"); + + UnmatchedME0Muon_ScatterPlot->Write(); UnmatchedME0Muon_ScatterPlot->Draw(); UnmatchedME0Muon_ScatterPlot->Print(histoFolder+"/UnmatchedME0Muon_ScatterPlot.png"); + + //gStyle->SetOptStat('oue'); c1->SetLogy(); UnmatchedME0Muon_Pt->Write(); UnmatchedME0Muon_Pt->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Pt.png"); + UnmatchedME0Muon_SmallBins_Pt->Write(); UnmatchedME0Muon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_SmallBins_Pt.png"); + UnmatchedME0Muon_VariableBins_Pt->Write(); UnmatchedME0Muon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_VariableBins_Pt.png"); + + Chi2UnmatchedME0Muon_Pt->Write(); Chi2UnmatchedME0Muon_Pt->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_Pt.png"); + Chi2UnmatchedME0Muon_SmallBins_Pt->Write(); Chi2UnmatchedME0Muon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_SmallBins_Pt.png"); + Chi2UnmatchedME0Muon_VariableBins_Pt->Write(); Chi2UnmatchedME0Muon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_VariableBins_Pt.png"); + + + Chi2MatchedME0Muon_Pt->Write(); Chi2MatchedME0Muon_Pt->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_Pt.png"); + Chi2MatchedME0Muon_SmallBins_Pt->Write(); Chi2MatchedME0Muon_SmallBins_Pt->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_SmallBins_Pt.png"); + Chi2MatchedME0Muon_VariableBins_Pt->Write(); Chi2MatchedME0Muon_VariableBins_Pt->Draw(); c1->Print(histoFolder+"/Chi2MatchedME0Muon_VariableBins_Pt.png"); + + UnmatchedME0Muon_Window_Pt->Write(); UnmatchedME0Muon_Window_Pt->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Window_Pt.png"); gStyle->SetOptStat(0); + + FailedTrack_Window_XPull->Write(); FailedTrack_Window_XPull->Draw(); c1->Print(histoFolder+"/FailedTrack_Window_XPull.png"); + FailedTrack_Window_YPull->Write(); FailedTrack_Window_YPull->Draw(); c1->Print(histoFolder+"/FailedTrack_Window_YPull.png"); + FailedTrack_Window_XDiff->Write(); FailedTrack_Window_XDiff->Draw(); c1->Print(histoFolder+"/FailedTrack_Window_XDiff.png"); + FailedTrack_Window_YDiff->Write(); FailedTrack_Window_YDiff->Draw(); c1->Print(histoFolder+"/FailedTrack_Window_YDiff.png"); + FailedTrack_Window_PhiDiff->Write(); FailedTrack_Window_PhiDiff->Draw(); c1->Print(histoFolder+"/FailedTrack_Window_PhiDiff.png"); + c1->SetLogy(0); TH1F *UnmatchedME0Muon_Cuts_Eta_PerEvent; - UnmatchedME0Muon_Cuts_Eta_PerEvent = new TH1F("UnmatchedME0Muon_Cuts_Eta_PerEvent" , "Muon |#eta|" , 4, 2.0, 2.8 ); + UnmatchedME0Muon_Cuts_Eta_PerEvent = new TH1F("UnmatchedME0Muon_Cuts_Eta_PerEvent" , "Muon |#eta|" , 8, 2.0, 2.8 ); //UnmatchedME0Muon_Cuts_Eta_PerEvent->Sumw2(); for (int i=1; i<=UnmatchedME0Muon_Cuts_Eta_PerEvent->GetNbinsX(); ++i){ UnmatchedME0Muon_Cuts_Eta_PerEvent->SetBinContent(i,(UnmatchedME0Muon_Cuts_Eta->GetBinContent(i))); @@ -1610,8 +2405,26 @@ void ME0MuonAnalyzer::endJob() UnmatchedME0Muon_Cuts_Eta_PerEvent->Write(); UnmatchedME0Muon_Cuts_Eta_PerEvent->Draw(); c1->Print(histoFolder+"/UnmatchedME0Muon_Cuts_Eta_PerEvent.png"); + + TH1F *Chi2UnmatchedME0Muon_Eta_PerEvent; + Chi2UnmatchedME0Muon_Eta_PerEvent = new TH1F("Chi2UnmatchedME0Muon_Eta_PerEvent" , "Muon |#eta|" , 8, 2.0, 2.8 ); + //Chi2UnmatchedME0Muon_Eta_PerEvent->Sumw2(); + for (int i=1; i<=Chi2UnmatchedME0Muon_Eta_PerEvent->GetNbinsX(); ++i){ + Chi2UnmatchedME0Muon_Eta_PerEvent->SetBinContent(i,(Chi2UnmatchedME0Muon_Eta->GetBinContent(i))); + } + Chi2UnmatchedME0Muon_Eta_PerEvent->Scale(1/Nevents); + + Chi2UnmatchedME0Muon_Eta_PerEvent->GetXaxis()->SetTitle("ME0Muon |#eta|"); + Chi2UnmatchedME0Muon_Eta_PerEvent->GetXaxis()->SetTitleSize(0.05); + + Chi2UnmatchedME0Muon_Eta_PerEvent->GetYaxis()->SetTitle("Average \# ME0Muons per event"); + Chi2UnmatchedME0Muon_Eta_PerEvent->GetYaxis()->SetTitleSize(0.05); + + Chi2UnmatchedME0Muon_Eta_PerEvent->Write(); Chi2UnmatchedME0Muon_Eta_PerEvent->Draw(); c1->Print(histoFolder+"/Chi2UnmatchedME0Muon_Eta_PerEvent.png"); + + TH1F *ME0Muon_Cuts_Eta_PerEvent; - ME0Muon_Cuts_Eta_PerEvent = new TH1F("ME0Muon_Cuts_Eta_PerEvent" , "Muon |#eta|" , 4, 2.0, 2.8 ); + ME0Muon_Cuts_Eta_PerEvent = new TH1F("ME0Muon_Cuts_Eta_PerEvent" , "Muon |#eta|" , 8, 2.0, 2.8 ); for (int i=1; i<=ME0Muon_Cuts_Eta_PerEvent->GetNbinsX(); ++i){ ME0Muon_Cuts_Eta_PerEvent->SetBinContent(i,(ME0Muon_Cuts_Eta->GetBinContent(i))); @@ -1637,6 +2450,21 @@ void ME0MuonAnalyzer::endJob() ClosestDelR_s->SetMarkerSize(3.0); ClosestDelR_s->Write(); ClosestDelR_s->Draw(); c1->Print(histoFolder+"/ClosestDelR_s.png"); + DelR_Window_Under5->Write(); DelR_Window_Under5->Draw(); c1->Print(histoFolder+"/DelR_Window_Under5.png"); + Pt_Window_Under5->Write(); Pt_Window_Under5->Draw(); c1->Print(histoFolder+"/Pt_Window_Under5.png"); + + DelR_Track_Window_Under5->Write(); DelR_Track_Window_Under5->Draw(); c1->Print(histoFolder+"/DelR_Track_Window_Under5.png"); + Pt_Track_Window_Under5->Write(); Pt_Track_Window_Under5->Draw(); c1->Print(histoFolder+"/Pt_Track_Window_Under5.png"); + c1->SetLogy(1); + Pt_Track_Window->Write(); Pt_Track_Window->Draw(); c1->Print(histoFolder+"/Pt_Track_Window.png"); + c1->SetLogy(0); + + DelR_Track_Window_Failed_Under5->Write(); DelR_Track_Window_Failed_Under5->Draw(); c1->Print(histoFolder+"/DelR_Track_Window_Failed_Under5.png"); + Pt_Track_Window_Failed_Under5->Write(); Pt_Track_Window_Failed_Under5->Draw(); c1->Print(histoFolder+"/Pt_Track_Window_Failed_Under5.png"); + c1->SetLogy(1); + Pt_Track_Window_Failed->Write(); Pt_Track_Window_Failed->Draw(); c1->Print(histoFolder+"/Pt_Track_Window_Failed.png"); + c1->SetLogy(0); + DelR_Segment_GenMuon->Write(); DelR_Segment_GenMuon->Draw(); c1->Print(histoFolder+"/DelR_Segment_GenMuon.png"); ClosestDelR_p->GetXaxis()->SetTitle("Gen Muon #eta"); @@ -1667,21 +2495,59 @@ void ME0MuonAnalyzer::endJob() FakeTracksPerAssociatedSegment_p->GetYaxis()->SetTitle("Average N_{Tracks} per segment"); FakeTracksPerAssociatedSegment_p->Write(); FakeTracksPerAssociatedSegment_p->Draw(); c1->Print(histoFolder+"/FakeTracksPerAssociatedSegment_p.png"); - GenMuon_Eta->Sumw2(); MatchedME0Muon_Eta->Sumw2(); Chi2MatchedME0Muon_Eta->Sumw2(); Chi2UnmatchedME0Muon_Eta->Sumw2(); - GenMuon_Pt->Sumw2(); MatchedME0Muon_Pt->Sumw2(); + PreMatch_TP_R->Write(); PreMatch_TP_R->Draw(); c1->Print(histoFolder+"/PreMatch_TP_R.png"); + PostMatch_TP_R->Write(); PostMatch_TP_R->Draw(); c1->Print(histoFolder+"/PostMatch_TP_R.png"); + PostMatch_BX0_TP_R->Write(); PostMatch_BX0_TP_R->Draw(); c1->Print(histoFolder+"/PostMatch_BX0_TP_R.png"); + + GenMuon_Eta->Sumw2(); MatchedME0Muon_Eta->Sumw2(); Chi2MatchedME0Muon_Eta->Sumw2(); Chi2UnmatchedME0Muon_Eta->Sumw2();TPMuon_Eta->Sumw2(); + GenMuon_Pt->Sumw2(); MatchedME0Muon_Pt->Sumw2(); MatchedME0Muon_SmallBins_Pt->Sumw2(); MatchedME0Muon_VariableBins_Pt->Sumw2(); + StandaloneMatchedME0Muon_Eta->Sumw2(); + StandaloneMatchedME0Muon_WideBinning_Eta->Sumw2(); + StandaloneMatchedME0Muon_WidestBinning_Eta->Sumw2(); + Track_Eta->Sumw2(); ME0Muon_Eta->Sumw2(); - Track_Pt->Sumw2(); ME0Muon_Pt->Sumw2(); + Track_Pt->Sumw2(); ME0Muon_Pt->Sumw2(); ME0Muon_SmallBins_Pt->Sumw2(); ME0Muon_VariableBins_Pt->Sumw2(); UnmatchedME0Muon_Eta->Sumw2(); UnmatchedME0Muon_Pt->Sumw2(); + UnmatchedME0Muon_SmallBins_Pt->Sumw2(); UnmatchedME0Muon_VariableBins_Pt->Sumw2(); UnmatchedME0Muon_Cuts_Eta->Sumw2(); ME0Muon_Cuts_Eta->Sumw2(); - ME0Muon_Cuts_Eta_5_10->Sumw2(); ME0Muon_Cuts_Eta_10_20->Sumw2(); ME0Muon_Cuts_Eta_20_40->Sumw2(); ME0Muon_Cuts_Eta_40->Sumw2(); - UnmatchedME0Muon_Cuts_Eta_5_10->Sumw2(); UnmatchedME0Muon_Cuts_Eta_10_20->Sumw2(); UnmatchedME0Muon_Cuts_Eta_20_40->Sumw2(); UnmatchedME0Muon_Cuts_Eta_40->Sumw2(); - GenMuon_Eta_5_10->Sumw2(); GenMuon_Eta_10_20->Sumw2(); GenMuon_Eta_20_40->Sumw2(); GenMuon_Eta_40->Sumw2(); - MatchedME0Muon_Eta_5_10->Sumw2(); MatchedME0Muon_Eta_10_20->Sumw2(); MatchedME0Muon_Eta_20_40->Sumw2(); MatchedME0Muon_Eta_40->Sumw2(); + ME0Muon_Cuts_Eta_5_10->Sumw2(); ME0Muon_Cuts_Eta_9_11->Sumw2(); ME0Muon_Cuts_Eta_10_50->Sumw2(); ME0Muon_Cuts_Eta_50_100->Sumw2(); ME0Muon_Cuts_Eta_100->Sumw2(); + UnmatchedME0Muon_Cuts_Eta_5_10->Sumw2(); UnmatchedME0Muon_Cuts_Eta_9_11->Sumw2(); UnmatchedME0Muon_Cuts_Eta_10_50->Sumw2(); UnmatchedME0Muon_Cuts_Eta_50_100->Sumw2(); UnmatchedME0Muon_Cuts_Eta_100->Sumw2(); + GenMuon_Eta_5_10->Sumw2(); GenMuon_Eta_9_11->Sumw2(); GenMuon_Eta_10_20->Sumw2(); GenMuon_Eta_20_40->Sumw2(); GenMuon_Eta_40->Sumw2(); + MatchedME0Muon_Eta_5_10->Sumw2(); MatchedME0Muon_Eta_9_11->Sumw2(); MatchedME0Muon_Eta_10_20->Sumw2(); MatchedME0Muon_Eta_20_40->Sumw2(); MatchedME0Muon_Eta_40->Sumw2(); + + Chi2MatchedME0Muon_Eta_5_10->Sumw2(); Chi2MatchedME0Muon_Eta_9_11->Sumw2(); Chi2MatchedME0Muon_Eta_10_20->Sumw2(); Chi2MatchedME0Muon_Eta_20_40->Sumw2(); Chi2MatchedME0Muon_Eta_40->Sumw2(); + TPMuon_Eta_5_10->Sumw2(); TPMuon_Eta_10_20->Sumw2(); TPMuon_Eta_20_40->Sumw2(); TPMuon_Eta_40->Sumw2(); + + + ME0Muon_Cuts_Eta_5_10->Write(); ME0Muon_Cuts_Eta_9_11->Write(); ME0Muon_Cuts_Eta_10_50->Write(); ME0Muon_Cuts_Eta_50_100->Write(); ME0Muon_Cuts_Eta_100->Write(); + UnmatchedME0Muon_Cuts_Eta_5_10->Write(); UnmatchedME0Muon_Cuts_Eta_9_11->Write(); UnmatchedME0Muon_Cuts_Eta_10_50->Write(); UnmatchedME0Muon_Cuts_Eta_50_100->Write(); UnmatchedME0Muon_Cuts_Eta_100->Write(); + GenMuon_Eta_5_10->Write(); GenMuon_Eta_9_11->Write(); GenMuon_Eta_10_20->Write(); GenMuon_Eta_20_40->Write(); GenMuon_Eta_40->Write(); + MatchedME0Muon_Eta_5_10->Write(); MatchedME0Muon_Eta_9_11->Write(); MatchedME0Muon_Eta_10_20->Write(); MatchedME0Muon_Eta_20_40->Write(); MatchedME0Muon_Eta_40->Write(); + + Chi2MatchedME0Muon_Eta_5_10->Write(); Chi2MatchedME0Muon_Eta_9_11->Write(); Chi2MatchedME0Muon_Eta_10_20->Write(); Chi2MatchedME0Muon_Eta_20_40->Write(); Chi2MatchedME0Muon_Eta_40->Write(); + TPMuon_Eta_5_10->Write(); TPMuon_Eta_10_20->Write(); TPMuon_Eta_20_40->Write(); TPMuon_Eta_40->Write(); + + UnmatchedME0Muon_Cuts_WideBinning_Eta->Sumw2(); + UnmatchedME0Muon_Cuts_WidestBinning_Eta->Sumw2(); + GenMuon_WideBinning_Eta->Sumw2(); + GenMuon_WidestBinning_Eta->Sumw2(); + TPMuon_WideBinning_Eta->Sumw2(); + TPMuon_WidestBinning_Eta->Sumw2(); + MatchedME0Muon_WideBinning_Eta->Sumw2(); + MatchedME0Muon_WidestBinning_Eta->Sumw2(); + Chi2MatchedME0Muon_WideBinning_Eta->Sumw2(); + Chi2MatchedME0Muon_WidestBinning_Eta->Sumw2(); + ME0Muon_Cuts_WideBinning_Eta->Sumw2(); + ME0Muon_Cuts_WidestBinning_Eta->Sumw2(); + Chi2UnmatchedME0Muon_WideBinning_Eta->Sumw2(); + Chi2UnmatchedME0Muon_WidestBinning_Eta->Sumw2(); + + //Captions/labels std::stringstream PtCutString; @@ -1742,6 +2608,96 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/MuonRecoEff_Eta.png"); + MuonRecoEff_WideBinning_Eta->Divide(MatchedME0Muon_WideBinning_Eta, GenMuon_WideBinning_Eta, 1, 1, "B"); + MuonRecoEff_WideBinning_Eta->GetXaxis()->SetTitle("Gen Muon |#eta|"); + MuonRecoEff_WideBinning_Eta->GetXaxis()->SetTitleSize(0.05); + MuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitle("ME0Muon Efficiency"); + MuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //MuonRecoEff_WideBinning_Eta->SetMinimum(MuonRecoEff_WideBinning_Eta->GetMinimum()-0.1); + MuonRecoEff_WideBinning_Eta->SetMinimum(0); + //MuonRecoEff_WideBinning_Eta->SetMaximum(MuonRecoEff_WideBinning_Eta->GetMaximum()+0.1); + MuonRecoEff_WideBinning_Eta->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + MuonRecoEff_WideBinning_Eta->Write(); MuonRecoEff_WideBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestMuonRecoEff_WideBinning_Eta.png"); + c1->Print(histoFolder+"/MuonRecoEff_WideBinning_Eta.png"); + + + MuonRecoEff_WidestBinning_Eta->Divide(MatchedME0Muon_WidestBinning_Eta, GenMuon_WidestBinning_Eta, 1, 1, "B"); + MuonRecoEff_WidestBinning_Eta->GetXaxis()->SetTitle("Gen Muon |#eta|"); + MuonRecoEff_WidestBinning_Eta->GetXaxis()->SetTitleSize(0.05); + MuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitle("ME0Muon Efficiency"); + MuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //MuonRecoEff_WidestBinning_Eta->SetMinimum(MuonRecoEff_WidestBinning_Eta->GetMinimum()-0.1); + MuonRecoEff_WidestBinning_Eta->SetMinimum(0); + //MuonRecoEff_WidestBinning_Eta->SetMaximum(MuonRecoEff_WidestBinning_Eta->GetMaximum()+0.1); + MuonRecoEff_WidestBinning_Eta->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + MuonRecoEff_WidestBinning_Eta->Write(); MuonRecoEff_WidestBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestMuonRecoEff_WidestBinning_Eta.png"); + c1->Print(histoFolder+"/MuonRecoEff_WidestBinning_Eta.png"); + + + StandaloneMuonRecoEff_Eta->Divide(StandaloneMatchedME0Muon_Eta, GenMuon_Eta, 1, 1, "B"); + StandaloneMuonRecoEff_Eta->GetXaxis()->SetTitle("Gen Muon |#eta|"); + StandaloneMuonRecoEff_Eta->GetXaxis()->SetTitleSize(0.05); + StandaloneMuonRecoEff_Eta->GetYaxis()->SetTitle("Standalone Muon Efficiency"); + StandaloneMuonRecoEff_Eta->GetYaxis()->SetTitleSize(0.05); + //StandaloneMuonRecoEff_Eta->SetMinimum(StandaloneMuonRecoEff_Eta->GetMinimum()-0.1); + StandaloneMuonRecoEff_Eta->SetMinimum(0); + //StandaloneMuonRecoEff_Eta->SetMaximum(StandaloneMuonRecoEff_Eta->GetMaximum()+0.1); + StandaloneMuonRecoEff_Eta->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + StandaloneMuonRecoEff_Eta->Write(); StandaloneMuonRecoEff_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestStandaloneMuonRecoEff_Eta.png"); + c1->Print(histoFolder+"/StandaloneMuonRecoEff_Eta.png"); + + + StandaloneMuonRecoEff_WideBinning_Eta->Divide(StandaloneMatchedME0Muon_WideBinning_Eta, GenMuon_WideBinning_Eta, 1, 1, "B"); + StandaloneMuonRecoEff_WideBinning_Eta->GetXaxis()->SetTitle("Gen Muon |#eta|"); + StandaloneMuonRecoEff_WideBinning_Eta->GetXaxis()->SetTitleSize(0.05); + StandaloneMuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitle("Standalone Muon Efficiency"); + StandaloneMuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //StandaloneMuonRecoEff_WideBinning_Eta->SetMinimum(StandaloneMuonRecoEff_WideBinning_Eta->GetMinimum()-0.1); + StandaloneMuonRecoEff_WideBinning_Eta->SetMinimum(0); + //StandaloneMuonRecoEff_WideBinning_Eta->SetMaximum(StandaloneMuonRecoEff_WideBinning_Eta->GetMaximum()+0.1); + StandaloneMuonRecoEff_WideBinning_Eta->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + StandaloneMuonRecoEff_WideBinning_Eta->Write(); StandaloneMuonRecoEff_WideBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestStandaloneMuonRecoEff_WideBinning_Eta.png"); + c1->Print(histoFolder+"/StandaloneMuonRecoEff_WideBinning_Eta.png"); + + + StandaloneMuonRecoEff_WidestBinning_Eta->Divide(StandaloneMatchedME0Muon_WidestBinning_Eta, GenMuon_WidestBinning_Eta, 1, 1, "B"); + StandaloneMuonRecoEff_WidestBinning_Eta->GetXaxis()->SetTitle("Gen Muon |#eta|"); + StandaloneMuonRecoEff_WidestBinning_Eta->GetXaxis()->SetTitleSize(0.05); + StandaloneMuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitle("Standalone Muon Efficiency"); + StandaloneMuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //StandaloneMuonRecoEff_WidestBinning_Eta->SetMinimum(StandaloneMuonRecoEff_WidestBinning_Eta->GetMinimum()-0.1); + StandaloneMuonRecoEff_WidestBinning_Eta->SetMinimum(0); + //StandaloneMuonRecoEff_WidestBinning_Eta->SetMaximum(StandaloneMuonRecoEff_WidestBinning_Eta->GetMaximum()+0.1); + StandaloneMuonRecoEff_WidestBinning_Eta->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + StandaloneMuonRecoEff_WidestBinning_Eta->Write(); StandaloneMuonRecoEff_WidestBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestStandaloneMuonRecoEff_WidestBinning_Eta.png"); + c1->Print(histoFolder+"/StandaloneMuonRecoEff_WidestBinning_Eta.png"); + + MuonRecoEff_Eta_5_10->Divide(MatchedME0Muon_Eta_5_10, GenMuon_Eta_5_10, 1, 1, "B"); MuonRecoEff_Eta_5_10->GetXaxis()->SetTitle("Gen Muon |#eta|"); MuonRecoEff_Eta_5_10->GetXaxis()->SetTitleSize(0.05); @@ -1759,6 +2715,23 @@ void ME0MuonAnalyzer::endJob() //c1->SaveAs("TestMuonRecoEff_Eta_5_10.png"); c1->Print(histoFolder+"/MuonRecoEff_Eta_5_10.png"); + MuonRecoEff_Eta_9_11->Divide(MatchedME0Muon_Eta_9_11, GenMuon_Eta_9_11, 1, 1, "B"); + MuonRecoEff_Eta_9_11->GetXaxis()->SetTitle("Gen Muon |#eta|"); + MuonRecoEff_Eta_9_11->GetXaxis()->SetTitleSize(0.05); + MuonRecoEff_Eta_9_11->GetYaxis()->SetTitle("ME0Muon Efficiency"); + MuonRecoEff_Eta_9_11->GetYaxis()->SetTitleSize(0.05); + //MuonRecoEff_Eta_9_11->SetMinimum(MuonRecoEff_Eta_9_11->GetMinimum()-0.1); + MuonRecoEff_Eta_9_11->SetMinimum(0); + //MuonRecoEff_Eta_9_11->SetMaximum(MuonRecoEff_Eta_9_11->GetMaximum()+0.1); + MuonRecoEff_Eta_9_11->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + MuonRecoEff_Eta_9_11->Write(); MuonRecoEff_Eta_9_11->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestMuonRecoEff_Eta_9_11.png"); + c1->Print(histoFolder+"/MuonRecoEff_Eta_9_11.png"); + MuonRecoEff_Eta_10_20->Divide(MatchedME0Muon_Eta_10_20, GenMuon_Eta_10_20, 1, 1, "B"); MuonRecoEff_Eta_10_20->GetXaxis()->SetTitle("Gen Muon |#eta|"); MuonRecoEff_Eta_10_20->GetXaxis()->SetTitleSize(0.05); @@ -1813,8 +2786,11 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/MuonRecoEff_Eta_40.png"); - Chi2MuonRecoEff_Eta->Divide(Chi2MatchedME0Muon_Eta, GenMuon_Eta, 1, 1, "B"); - std::cout<<"GenMuon_Eta = "<Integral()<Divide(Chi2MatchedME0Muon_Eta, TPMuon_Eta, 1, 1, "B"); + std::cout<<"TPMuon_Eta = "<Integral()<GetXaxis()->SetTitleSize(0.05); @@ -1836,6 +2812,146 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta.png"); + Chi2MuonRecoEff_WideBinning_Eta->Divide(Chi2MatchedME0Muon_WideBinning_Eta, TPMuon_WideBinning_Eta, 1, 1, "B"); + std::cout<<"TPMuon_WideBinning_Eta = "<Integral()<GetXaxis()->SetTitleSize(0.05); + + Chi2MuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_WideBinning_Eta->GetYaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_WideBinning_Eta->SetMinimum(0); + Chi2MuonRecoEff_WideBinning_Eta->SetMaximum(1.2); + + Chi2MuonRecoEff_WideBinning_Eta->Write(); Chi2MuonRecoEff_WideBinning_Eta->Draw(); + + txt->DrawLatex(0.15,0.2,pcstr); + + + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //latex->SetTextAlign(align_); + latex->DrawLatex(0.4, 0.85, cmsText); + + c1->Print(histoFolder+"/Chi2MuonRecoEff_WideBinning_Eta.png"); + + Chi2MuonRecoEff_WidestBinning_Eta->Divide(Chi2MatchedME0Muon_WidestBinning_Eta, TPMuon_WidestBinning_Eta, 1, 1, "B"); + std::cout<<"TPMuon_WidestBinning_Eta = "<Integral()<GetXaxis()->SetTitleSize(0.05); + + Chi2MuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_WidestBinning_Eta->GetYaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_WidestBinning_Eta->SetMinimum(0); + Chi2MuonRecoEff_WidestBinning_Eta->SetMaximum(1.2); + + Chi2MuonRecoEff_WidestBinning_Eta->Write(); Chi2MuonRecoEff_WidestBinning_Eta->Draw(); + + txt->DrawLatex(0.15,0.2,pcstr); + + + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //latex->SetTextAlign(align_); + latex->DrawLatex(0.4, 0.85, cmsText); + + c1->Print(histoFolder+"/Chi2MuonRecoEff_WidestBinning_Eta.png"); + + std::cout<<"Here0"<Divide(Chi2MatchedME0Muon_Eta_5_10, GenMuon_Eta_5_10, 1, 1, "B"); + std::cout<<"Here0"<GetXaxis()->SetTitle("Gen Muon |#eta|"); + Chi2MuonRecoEff_Eta_5_10->GetXaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_Eta_5_10->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_Eta_5_10->GetYaxis()->SetTitleSize(0.05); + //Chi2MuonRecoEff_Eta_5_10->SetMinimum(Chi2MuonRecoEff_Eta_5_10->GetMinimum()-0.1); + Chi2MuonRecoEff_Eta_5_10->SetMinimum(0); + //Chi2MuonRecoEff_Eta_5_10->SetMaximum(Chi2MuonRecoEff_Eta_5_10->GetMaximum()+0.1); + Chi2MuonRecoEff_Eta_5_10->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + Chi2MuonRecoEff_Eta_5_10->Write(); Chi2MuonRecoEff_Eta_5_10->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestChi2MuonRecoEff_Eta_5_10.png"); + c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta_5_10.png"); + + + Chi2MuonRecoEff_Eta_9_11->Divide(Chi2MatchedME0Muon_Eta_9_11, GenMuon_Eta_9_11, 1, 1, "B"); + std::cout<<"Here0"<GetXaxis()->SetTitle("Gen Muon |#eta|"); + Chi2MuonRecoEff_Eta_9_11->GetXaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_Eta_9_11->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_Eta_9_11->GetYaxis()->SetTitleSize(0.05); + //Chi2MuonRecoEff_Eta_9_11->SetMinimum(Chi2MuonRecoEff_Eta_9_11->GetMinimum()-0.1); + Chi2MuonRecoEff_Eta_9_11->SetMinimum(0); + //Chi2MuonRecoEff_Eta_9_11->SetMaximum(Chi2MuonRecoEff_Eta_9_11->GetMaximum()+0.1); + Chi2MuonRecoEff_Eta_9_11->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + Chi2MuonRecoEff_Eta_9_11->Write(); Chi2MuonRecoEff_Eta_9_11->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestChi2MuonRecoEff_Eta_9_11.png"); + c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta_9_11.png"); + + std::cout<<"Here"<Divide(Chi2MatchedME0Muon_Eta_10_20, GenMuon_Eta_10_20, 1, 1, "B"); + Chi2MuonRecoEff_Eta_10_20->GetXaxis()->SetTitle("Gen Muon |#eta|"); + Chi2MuonRecoEff_Eta_10_20->GetXaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_Eta_10_20->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_Eta_10_20->GetYaxis()->SetTitleSize(0.05); + //Chi2MuonRecoEff_Eta_10_20->SetMinimum(Chi2MuonRecoEff_Eta_10_20->GetMinimum()-0.1); + Chi2MuonRecoEff_Eta_10_20->SetMinimum(0); + //Chi2MuonRecoEff_Eta_10_20->SetMaximum(Chi2MuonRecoEff_Eta_10_20->GetMaximum()+0.1); + Chi2MuonRecoEff_Eta_10_20->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + Chi2MuonRecoEff_Eta_10_20->Write(); Chi2MuonRecoEff_Eta_10_20->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestChi2MuonRecoEff_Eta_10_20.png"); + c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta_10_20.png"); + + + Chi2MuonRecoEff_Eta_20_40->Divide(Chi2MatchedME0Muon_Eta_20_40, GenMuon_Eta_20_40, 1, 1, "B"); + Chi2MuonRecoEff_Eta_20_40->GetXaxis()->SetTitle("Gen Muon |#eta|"); + Chi2MuonRecoEff_Eta_20_40->GetXaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_Eta_20_40->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_Eta_20_40->GetYaxis()->SetTitleSize(0.05); + //Chi2MuonRecoEff_Eta_20_40->SetMinimum(Chi2MuonRecoEff_Eta_20_40->GetMinimum()-0.1); + Chi2MuonRecoEff_Eta_20_40->SetMinimum(0); + //Chi2MuonRecoEff_Eta_20_40->SetMaximum(Chi2MuonRecoEff_Eta_20_40->GetMaximum()+0.1); + Chi2MuonRecoEff_Eta_20_40->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + Chi2MuonRecoEff_Eta_20_40->Write(); Chi2MuonRecoEff_Eta_20_40->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestChi2MuonRecoEff_Eta_20_40.png"); + c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta_20_40.png"); + + + Chi2MuonRecoEff_Eta_40->Divide(Chi2MatchedME0Muon_Eta_40, GenMuon_Eta_40, 1, 1, "B"); + Chi2MuonRecoEff_Eta_40->GetXaxis()->SetTitle("Gen Muon |#eta|"); + Chi2MuonRecoEff_Eta_40->GetXaxis()->SetTitleSize(0.05); + Chi2MuonRecoEff_Eta_40->GetYaxis()->SetTitle("ME0Muon Efficiency"); + Chi2MuonRecoEff_Eta_40->GetYaxis()->SetTitleSize(0.05); + //Chi2MuonRecoEff_Eta_40->SetMinimum(Chi2MuonRecoEff_Eta_40->GetMinimum()-0.1); + Chi2MuonRecoEff_Eta_40->SetMinimum(0); + //Chi2MuonRecoEff_Eta_40->SetMaximum(Chi2MuonRecoEff_Eta_40->GetMaximum()+0.1); + Chi2MuonRecoEff_Eta_40->SetMaximum(1.2); + //CMS_lumi( c1, 7, 11 ); + Chi2MuonRecoEff_Eta_40->Write(); Chi2MuonRecoEff_Eta_40->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + + //c1->SaveAs("TestChi2MuonRecoEff_Eta_40.png"); + c1->Print(histoFolder+"/Chi2MuonRecoEff_Eta_40.png"); + std::cout<<" MuonRecoEff_Eta values:"<Sumw2(); for (int i=1; i<=MuonRecoEff_Eta->GetNbinsX(); ++i){ @@ -1872,6 +2988,42 @@ void ME0MuonAnalyzer::endJob() //c1->SaveAs('TestFakeRate_Eta.png'); c1->Print(histoFolder+"/FakeRate_Eta.png"); + FakeRate_WideBinning_Eta->Divide(UnmatchedME0Muon_Cuts_WideBinning_Eta, ME0Muon_Cuts_WideBinning_Eta, 1, 1, "B"); + FakeRate_WideBinning_Eta->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_WideBinning_Eta->GetXaxis()->SetTitleSize(0.05); + FakeRate_WideBinning_Eta->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_WideBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //FakeRate_WideBinning_Eta->SetMinimum(FakeRate_WideBinning_Eta->GetMinimum()-0.1); + FakeRate_WideBinning_Eta->SetMinimum(0); + //FakeRate_WideBinning_Eta->SetMaximum(FakeRate_WideBinning_Eta->GetMaximum()+0.1); + FakeRate_WideBinning_Eta->SetMaximum(1.2); + FakeRate_WideBinning_Eta->Write(); FakeRate_WideBinning_Eta->Draw(); + + txt->DrawLatex(0.15,0.4,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //c1->SaveAs('TestFakeRate_WideBinning_Eta.png'); + c1->Print(histoFolder+"/FakeRate_WideBinning_Eta.png"); + + + FakeRate_WidestBinning_Eta->Divide(UnmatchedME0Muon_Cuts_WidestBinning_Eta, ME0Muon_Cuts_WidestBinning_Eta, 1, 1, "B"); + FakeRate_WidestBinning_Eta->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_WidestBinning_Eta->GetXaxis()->SetTitleSize(0.05); + FakeRate_WidestBinning_Eta->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_WidestBinning_Eta->GetYaxis()->SetTitleSize(0.05); + //FakeRate_WidestBinning_Eta->SetMinimum(FakeRate_WidestBinning_Eta->GetMinimum()-0.1); + FakeRate_WidestBinning_Eta->SetMinimum(0); + //FakeRate_WidestBinning_Eta->SetMaximum(FakeRate_WidestBinning_Eta->GetMaximum()+0.1); + FakeRate_WidestBinning_Eta->SetMaximum(1.2); + FakeRate_WidestBinning_Eta->Write(); FakeRate_WidestBinning_Eta->Draw(); + + txt->DrawLatex(0.15,0.4,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //c1->SaveAs('TestFakeRate_WidestBinning_Eta.png'); + c1->Print(histoFolder+"/FakeRate_WidestBinning_Eta.png"); FakeRate_Eta_5_10->Divide(UnmatchedME0Muon_Cuts_Eta_5_10, ME0Muon_Cuts_Eta_5_10, 1, 1, "B"); @@ -1893,64 +3045,83 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/FakeRate_Eta_5_10.png"); + FakeRate_Eta_9_11->Divide(UnmatchedME0Muon_Cuts_Eta_9_11, ME0Muon_Cuts_Eta_9_11, 1, 1, "B"); + FakeRate_Eta_9_11->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_Eta_9_11->GetXaxis()->SetTitleSize(0.05); + FakeRate_Eta_9_11->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_Eta_9_11->GetYaxis()->SetTitleSize(0.05); + //FakeRate_Eta_9_11->SetMinimum(FakeRate_Eta_9_11->GetMinimum()-0.1); + FakeRate_Eta_9_11->SetMinimum(0); + //FakeRate_Eta_9_11->SetMaximum(FakeRate_Eta_9_11->GetMaximum()+0.1); + FakeRate_Eta_9_11->SetMaximum(1.2); + FakeRate_Eta_9_11->Write(); FakeRate_Eta_9_11->Draw(); + + txt->DrawLatex(0.15,0.4,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //c1->SaveAs('TestFakeRate_Eta_9_11.png'); + c1->Print(histoFolder+"/FakeRate_Eta_9_11.png"); - FakeRate_Eta_10_20->Divide(UnmatchedME0Muon_Cuts_Eta_10_20, ME0Muon_Cuts_Eta_10_20, 1, 1, "B"); - FakeRate_Eta_10_20->GetXaxis()->SetTitle("Reconstructed track |#eta|"); - FakeRate_Eta_10_20->GetXaxis()->SetTitleSize(0.05); - FakeRate_Eta_10_20->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); - FakeRate_Eta_10_20->GetYaxis()->SetTitleSize(0.05); - //FakeRate_Eta_10_20->SetMinimum(FakeRate_Eta_10_20->GetMinimum()-0.1); - FakeRate_Eta_10_20->SetMinimum(0); - //FakeRate_Eta_10_20->SetMaximum(FakeRate_Eta_10_20->GetMaximum()+0.1); - FakeRate_Eta_10_20->SetMaximum(1.2); - FakeRate_Eta_10_20->Write(); FakeRate_Eta_10_20->Draw(); + + + FakeRate_Eta_10_50->Divide(UnmatchedME0Muon_Cuts_Eta_10_50, ME0Muon_Cuts_Eta_10_50, 1, 1, "B"); + FakeRate_Eta_10_50->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_Eta_10_50->GetXaxis()->SetTitleSize(0.05); + FakeRate_Eta_10_50->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_Eta_10_50->GetYaxis()->SetTitleSize(0.05); + //FakeRate_Eta_10_50->SetMinimum(FakeRate_Eta_10_50->GetMinimum()-0.1); + FakeRate_Eta_10_50->SetMinimum(0); + //FakeRate_Eta_10_50->SetMaximum(FakeRate_Eta_10_50->GetMaximum()+0.1); + FakeRate_Eta_10_50->SetMaximum(1.2); + FakeRate_Eta_10_50->Write(); FakeRate_Eta_10_50->Draw(); txt->DrawLatex(0.15,0.4,pcstr); latex->DrawLatex(0.4, 0.85, cmsText); latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); - //c1->SaveAs('TestFakeRate_Eta_10_20.png'); - c1->Print(histoFolder+"/FakeRate_Eta_10_20.png"); + //c1->SaveAs('TestFakeRate_Eta_10_50.png'); + c1->Print(histoFolder+"/FakeRate_Eta_10_50.png"); - FakeRate_Eta_20_40->Divide(UnmatchedME0Muon_Cuts_Eta_20_40, ME0Muon_Cuts_Eta_20_40, 1, 1, "B"); - FakeRate_Eta_20_40->GetXaxis()->SetTitle("Reconstructed track |#eta|"); - FakeRate_Eta_20_40->GetXaxis()->SetTitleSize(0.05); - FakeRate_Eta_20_40->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); - FakeRate_Eta_20_40->GetYaxis()->SetTitleSize(0.05); - //FakeRate_Eta_20_40->SetMinimum(FakeRate_Eta_20_40->GetMinimum()-0.1); - FakeRate_Eta_20_40->SetMinimum(0); - //FakeRate_Eta_20_40->SetMaximum(FakeRate_Eta_20_40->GetMaximum()+0.1); - FakeRate_Eta_20_40->SetMaximum(1.2); - FakeRate_Eta_20_40->Write(); FakeRate_Eta_20_40->Draw(); + FakeRate_Eta_50_100->Divide(UnmatchedME0Muon_Cuts_Eta_50_100, ME0Muon_Cuts_Eta_50_100, 1, 1, "B"); + FakeRate_Eta_50_100->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_Eta_50_100->GetXaxis()->SetTitleSize(0.05); + FakeRate_Eta_50_100->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_Eta_50_100->GetYaxis()->SetTitleSize(0.05); + //FakeRate_Eta_50_100->SetMinimum(FakeRate_Eta_50_100->GetMinimum()-0.1); + FakeRate_Eta_50_100->SetMinimum(0); + //FakeRate_Eta_50_100->SetMaximum(FakeRate_Eta_50_100->GetMaximum()+0.1); + FakeRate_Eta_50_100->SetMaximum(1.2); + FakeRate_Eta_50_100->Write(); FakeRate_Eta_50_100->Draw(); txt->DrawLatex(0.15,0.4,pcstr); latex->DrawLatex(0.4, 0.85, cmsText); latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); - //c1->SaveAs('TestFakeRate_Eta_20_40.png'); - c1->Print(histoFolder+"/FakeRate_Eta_20_40.png"); + //c1->SaveAs('TestFakeRate_Eta_50_100.png'); + c1->Print(histoFolder+"/FakeRate_Eta_50_100.png"); - FakeRate_Eta_40->Divide(UnmatchedME0Muon_Cuts_Eta_40, ME0Muon_Cuts_Eta_40, 1, 1, "B"); - FakeRate_Eta_40->GetXaxis()->SetTitle("Reconstructed track |#eta|"); - FakeRate_Eta_40->GetXaxis()->SetTitleSize(0.05); - FakeRate_Eta_40->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); - FakeRate_Eta_40->GetYaxis()->SetTitleSize(0.05); - //FakeRate_Eta_40->SetMinimum(FakeRate_Eta_40->GetMinimum()-0.1); - FakeRate_Eta_40->SetMinimum(0); - //FakeRate_Eta_40->SetMaximum(FakeRate_Eta_40->GetMaximum()+0.1); - FakeRate_Eta_40->SetMaximum(1.2); - FakeRate_Eta_40->Write(); FakeRate_Eta_40->Draw(); + FakeRate_Eta_100->Divide(UnmatchedME0Muon_Cuts_Eta_100, ME0Muon_Cuts_Eta_100, 1, 1, "B"); + FakeRate_Eta_100->GetXaxis()->SetTitle("Reconstructed track |#eta|"); + FakeRate_Eta_100->GetXaxis()->SetTitleSize(0.05); + FakeRate_Eta_100->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + FakeRate_Eta_100->GetYaxis()->SetTitleSize(0.05); + //FakeRate_Eta_100->SetMinimum(FakeRate_Eta_100->GetMinimum()-0.1); + FakeRate_Eta_100->SetMinimum(0); + //FakeRate_Eta_100->SetMaximum(FakeRate_Eta_100->GetMaximum()+0.1); + FakeRate_Eta_100->SetMaximum(1.2); + FakeRate_Eta_100->Write(); FakeRate_Eta_100->Draw(); txt->DrawLatex(0.15,0.4,pcstr); latex->DrawLatex(0.4, 0.85, cmsText); latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); - //c1->SaveAs('TestFakeRate_Eta_40.png'); - c1->Print(histoFolder+"/FakeRate_Eta_40.png"); + //c1->SaveAs('TestFakeRate_Eta_100.png'); + c1->Print(histoFolder+"/FakeRate_Eta_100.png"); @@ -1974,6 +3145,50 @@ void ME0MuonAnalyzer::endJob() //c1->SaveAs('TestChi2FakeRate_Eta.png'); c1->Print(histoFolder+"/Chi2FakeRate_Eta.png"); + + + Chi2FakeRate_WideBinning_Eta->Divide(Chi2UnmatchedME0Muon_WideBinning_Eta, ME0Muon_Cuts_WideBinning_Eta, 1, 1, "B"); + //std::cout<<"Chi2UnmatchedME0Muon_WideBinning_Eta = "<Integral()<GetXaxis()->SetTitleSize(0.05); + Chi2FakeRate_WideBinning_Eta->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + Chi2FakeRate_WideBinning_Eta->GetYaxis()->SetTitleSize(0.05); + Chi2FakeRate_WideBinning_Eta->SetMinimum(0); + Chi2FakeRate_WideBinning_Eta->SetMaximum(1.2); + Chi2FakeRate_WideBinning_Eta->Write(); Chi2FakeRate_WideBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //c1->SaveAs('TestChi2FakeRate_WideBinning_Eta.png'); + c1->Print(histoFolder+"/Chi2FakeRate_WideBinning_Eta.png"); + + + Chi2FakeRate_WidestBinning_Eta->Divide(Chi2UnmatchedME0Muon_WidestBinning_Eta, ME0Muon_Cuts_WidestBinning_Eta, 1, 1, "B"); + //std::cout<<"Chi2UnmatchedME0Muon_WidestBinning_Eta = "<Integral()<GetXaxis()->SetTitleSize(0.05); + Chi2FakeRate_WidestBinning_Eta->GetYaxis()->SetTitle("ME0 Muon Fake Rate"); + Chi2FakeRate_WidestBinning_Eta->GetYaxis()->SetTitleSize(0.05); + Chi2FakeRate_WidestBinning_Eta->SetMinimum(0); + Chi2FakeRate_WidestBinning_Eta->SetMaximum(1.2); + Chi2FakeRate_WidestBinning_Eta->Write(); Chi2FakeRate_WidestBinning_Eta->Draw(); + txt->DrawLatex(0.15,0.2,pcstr); + latex->DrawLatex(0.4, 0.85, cmsText); + latex1->DrawLatex(1-r,1-t+lumiTextOffset*t,lumiText); + + //c1->SaveAs('TestChi2FakeRate_WidestBinning_Eta.png'); + c1->Print(histoFolder+"/Chi2FakeRate_WidestBinning_Eta.png"); + + //Fake Rate per Event: FakeRate_Eta_PerEvent->Divide(UnmatchedME0Muon_Cuts_Eta_PerEvent, ME0Muon_Cuts_Eta_PerEvent, 1, 1, "B"); @@ -2036,16 +3251,22 @@ void ME0MuonAnalyzer::endJob() PtDiff_s->GetYaxis()->SetTitle("|(pt track-ptgen)/ptgen|"); PtDiff_s->Write(); PtDiff_s->Draw(); c1->Print(histoFolder+"/PtDiff_s.png"); + PtDiff_s_5_10->Write(); + PtDiff_s_10_20->Write(); + PtDiff_s_20_40->Write(); + PtDiff_s_40->Write(); for(int i=1; i<=PtDiff_p->GetNbinsX(); ++i) { PtDiff_rms->SetBinContent(i, PtDiff_p->GetBinError(i)); std::cout<<"Pt_rms = "<GetBinError(i)<GetNbinsX(); ++i) { - + + TH1D *test; + std::cout<<"Total integral is "<Integral()<GetNbinsX(); ++i) { + + std::stringstream tempstore; tempstore<ProjectionY("test",i,i,""); + std::cout<<"Bin = "<GetBinContent(i)<Fit(gaus_narrow,"R"); @@ -2068,13 +3292,17 @@ void ME0MuonAnalyzer::endJob() // std::cout<Fit(gaus_wide,"R"); + std::cout<<"Getting values"<GetParameter(2); Double_t e_w2 = gaus_wide->GetParError(2); + std::cout<<"Got values"<SetBinContent(i, n2); // PtDiff_gaus_narrow->SetBinError(i, e_n2); PtDiff_gaus_wide->SetBinContent(i, w2); @@ -2084,13 +3312,23 @@ void ME0MuonAnalyzer::endJob() TString FileName = "Bin"+thistemp+"Fit.png"; c1->Print(histoFolder+"/"+FileName); - test->Draw(); + //test->Draw(); + //delete test; + + //continue; + delete test; + test= new TH1D("test" , "pt resolution" , 200, -1.0, 1.0 ); + test->Draw(); // Redoing for pt 5 to 10 + std::cout<<"About to project"<ProjectionY("test",i,i,""); + std::cout<<"About to check, "<Integral()<Integral() < 1.0) continue; + std::cout<<"Running the 5-10 fit"<Fit(gaus_5_10,"R"); @@ -2105,6 +3343,8 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/"+FileName); delete test; + test= new TH1D("test" , "pt resolution" , 200, -1.0, 1.0 ); + test->Draw(); // Redoing for pt 10 to 20 PtDiff_s_10_20->ProjectionY("test",i,i,""); if (test->Integral() < 1.0) continue; @@ -2123,6 +3363,9 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/"+FileName); delete test; + + test= new TH1D("test" , "pt resolution" , 200, -1.0, 1.0 ); + test->Draw(); // Redoing for pt 20 to 40 PtDiff_s_20_40->ProjectionY("test",i,i,""); if (test->Integral() < 1.0) continue; @@ -2141,6 +3384,9 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/"+FileName); delete test; + + test= new TH1D("test" , "pt resolution" , 200, -1.0, 1.0 ); + test->Draw(); // Redoing for pt 40+ PtDiff_s_40->ProjectionY("test",i,i,""); if (test->Integral() < 1.0) continue; @@ -2159,7 +3405,28 @@ void ME0MuonAnalyzer::endJob() c1->Print(histoFolder+"/"+FileName); delete test; - //test->Clear(); + + test= new TH1D("test" , "pt resolution" , 200, -1.0, 1.0 ); + test->Draw(); + // Redoing for pt 40+ + StandalonePtDiff_s->ProjectionY("test",i,i,""); + if (test->Integral() < 1.0) continue; + + TF1 *Standalonegaus = new TF1("Standalonegaus","gaus",-.2,.2); + test->Fit(Standalonegaus,"R"); + + w2 = gaus_40->GetParameter(2); + e_w2 = gaus_40->GetParError(2); + + StandalonePtDiff_gaus->SetBinContent(i, w2); + StandalonePtDiff_gaus->SetBinError(i, e_w2); + + test->Draw(); + FileName = "Bin"+thistemp+"StandaloneFit.png"; + c1->Print(histoFolder+"/"+FileName); + + delete test; + //test->Clear(); } // PtDiff_gaus_narrow->SetMarkerStyle(22); // PtDiff_gaus_narrow->SetMarkerSize(1.2); @@ -2192,7 +3459,7 @@ void ME0MuonAnalyzer::endJob() //PtDiff_gaus_5_10->Draw("PL"); PtDiff_gaus_5_10->GetXaxis()->SetTitle("Gen Muon |#eta|"); - PtDiff_gaus_5_10->GetYaxis()->SetTitle("Gaussian width of (pt track-ptgen)/ptgen"); + PtDiff_gaus_5_10->GetYaxis()->SetTitle("Gaussian width of (qoverpt track-qopverptp gen)/qoverpt gen"); PtDiff_gaus_5_10->Write(); PtDiff_gaus_5_10->Draw("PE"); c1->Print(histoFolder+"/PtDiff_gaus_5_10.png"); PtDiff_gaus_10_20->SetMarkerStyle(22); @@ -2203,7 +3470,7 @@ void ME0MuonAnalyzer::endJob() //PtDiff_gaus_10_20->Draw("PL"); PtDiff_gaus_10_20->GetXaxis()->SetTitle("Gen Muon |#eta|"); - PtDiff_gaus_10_20->GetYaxis()->SetTitle("Gaussian width of (pt track-ptgen)/ptgen"); + PtDiff_gaus_10_20->GetYaxis()->SetTitle("Gaussian width of (qoverpt track-qopverptp gen)/qoverpt gen"); PtDiff_gaus_10_20->Write(); PtDiff_gaus_10_20->Draw("PE"); c1->Print(histoFolder+"/PtDiff_gaus_10_20.png"); PtDiff_gaus_20_40->SetMarkerStyle(22); @@ -2214,7 +3481,7 @@ void ME0MuonAnalyzer::endJob() //PtDiff_gaus_20_40->Draw("PL"); PtDiff_gaus_20_40->GetXaxis()->SetTitle("Gen Muon |#eta|"); - PtDiff_gaus_20_40->GetYaxis()->SetTitle("Gaussian width of (pt track-ptgen)/ptgen"); + PtDiff_gaus_20_40->GetYaxis()->SetTitle("Gaussian width of (qoverpt track-qopverptp gen)/qoverpt gen"); PtDiff_gaus_20_40->Write(); PtDiff_gaus_20_40->Draw("PE"); c1->Print(histoFolder+"/PtDiff_gaus_20_40.png"); PtDiff_gaus_40->SetMarkerStyle(22); @@ -2225,9 +3492,20 @@ void ME0MuonAnalyzer::endJob() //PtDiff_gaus_40->Draw("PL"); PtDiff_gaus_40->GetXaxis()->SetTitle("Gen Muon |#eta|"); - PtDiff_gaus_40->GetYaxis()->SetTitle("Gaussian width of (pt track-ptgen)/ptgen"); + PtDiff_gaus_40->GetYaxis()->SetTitle("Gaussian width of (qoverpt track-qopverptp gen)/qoverpt gen"); PtDiff_gaus_40->Write(); PtDiff_gaus_40->Draw("PE"); c1->Print(histoFolder+"/PtDiff_gaus_40.png"); + StandalonePtDiff_gaus->SetMarkerStyle(22); + StandalonePtDiff_gaus->SetMarkerSize(1.2); + StandalonePtDiff_gaus->SetMarkerColor(kBlue); + //StandalonePtDiff_gaus->SetLineColor(kRed); + + //StandalonePtDiff_gaus->Draw("PL"); + + StandalonePtDiff_gaus->GetXaxis()->SetTitle("Gen Muon |#eta|"); + StandalonePtDiff_gaus->GetYaxis()->SetTitle("Gaussian width of (pt track-ptgen)/ptgen"); + StandalonePtDiff_gaus->Write(); StandalonePtDiff_gaus->Draw("PE"); c1->Print(histoFolder+"/StandalonePtDiff_gaus.png"); + PtDiff_p->SetMarkerStyle(1); PtDiff_p->SetMarkerSize(3.0); @@ -2282,6 +3560,11 @@ void ME0MuonAnalyzer::endJob() ofstream logout; logout.open (histoFolder+"/Log.txt"); + logout<<"Chi 2 Efficiencies and errors:\n"; + for (int i=1; i<=Chi2MuonRecoEff_Eta->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + logout<<"Efficiencies and errors:\n"; for (int i=1; i<=MuonRecoEff_Eta->GetNbinsX(); ++i){ logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; @@ -2303,6 +3586,18 @@ void ME0MuonAnalyzer::endJob() logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } + + logout<<"Efficiencies and errors 9_11:\n"; + for (int i=1; i<=MuonRecoEff_Eta_9_11->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + + + logout<<"Chi 2 Efficiencies and errors 5_10:\n"; + for (int i=1; i<=Chi2MuonRecoEff_Eta_5_10->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + logout<<"Fake Rate 5_10:\n"; for (int i=1; i<=FakeRate_Eta_5_10->GetNbinsX(); ++i){ logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; @@ -2319,9 +3614,15 @@ void ME0MuonAnalyzer::endJob() logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } - logout<<"Fake Rate 10_20:\n"; - for (int i=1; i<=FakeRate_Eta_10_20->GetNbinsX(); ++i){ - logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + logout<<"Chi 2 Efficiencies and errors 10_20:\n"; + for (int i=1; i<=Chi2MuonRecoEff_Eta_10_20->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + + + logout<<"Fake Rate 10_50:\n"; + for (int i=1; i<=FakeRate_Eta_10_50->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } logout<<"Resolution vs eta 10_20:\n"; @@ -2335,11 +3636,17 @@ void ME0MuonAnalyzer::endJob() logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } - logout<<"Fake Rate 20_40:\n"; - for (int i=1; i<=FakeRate_Eta_20_40->GetNbinsX(); ++i){ - logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + + logout<<"Chi 2 Efficiencies and errors 20_40:\n"; + for (int i=1; i<=Chi2MuonRecoEff_Eta_20_40->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } + logout<<"Fake Rate 50_100:\n"; + for (int i=1; i<=FakeRate_Eta_50_100->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + logout<<"Resolution vs eta 20_40:\n"; for (int i=1; i<=PtDiff_gaus_20_40->GetNbinsX(); ++i){ logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; @@ -2351,9 +3658,15 @@ void ME0MuonAnalyzer::endJob() logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } + + logout<<"Chi 2 Efficiencies and errors 40:\n"; + for (int i=1; i<=Chi2MuonRecoEff_Eta_40->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + } + logout<<"Fake Rate 40:\n"; - for (int i=1; i<=FakeRate_Eta_40->GetNbinsX(); ++i){ - logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; + for (int i=1; i<=FakeRate_Eta_100->GetNbinsX(); ++i){ + logout<GetBinContent(i)<<","<GetBinError(i)<<"\n"; } logout<<"Resolution vs eta 40:\n"; diff --git a/RecoMuon/MuonIdentification/test/testME0MuonAnalyzer_cfg.py b/RecoMuon/MuonIdentification/test/testME0MuonAnalyzer_cfg.py index c59229c4a1a09..18fef24c1fd66 100644 --- a/RecoMuon/MuonIdentification/test/testME0MuonAnalyzer_cfg.py +++ b/RecoMuon/MuonIdentification/test/testME0MuonAnalyzer_cfg.py @@ -36,7 +36,7 @@ process.source = cms.Source("PoolSource", - fileNames = cms.untracked.vstring('file:///somewhere/simevent.root') ##/somewhere/simevent.root" } + fileNames = cms.untracked.vstring('file:/store/group/upgrade/muon/ME0GlobalReco/ME0MuonReRun_DY_SLHC23patch1/M-20_TuneZ2star_14TeV_6_2_0_SLHC23patch1_2023/output_9_1_hVY.root') ) diff --git a/Validation/RecoMuon/plugins/ME0MuonTrackCollProducer.cc b/Validation/RecoMuon/plugins/ME0MuonTrackCollProducer.cc index c5192c415e37b..d43b5c01e31bb 100644 --- a/Validation/RecoMuon/plugins/ME0MuonTrackCollProducer.cc +++ b/Validation/RecoMuon/plugins/ME0MuonTrackCollProducer.cc @@ -55,12 +55,12 @@ ME0MuonTrackCollProducer::ME0MuonTrackCollProducer(const edm::ParameterSet& parset) : - muonsTag(parset.getParameter< edm::InputTag >("muonsTag")), - vxtTag(parset.getParameter< edm::InputTag >("vxtTag")), - useIPxy(parset.getUntrackedParameter< bool >("useIPxy", true)), - useIPz(parset.getUntrackedParameter< bool >("useIPz", true)), - selectionTags(parset.getParameter< std::vector >("selectionTags")), - trackType(parset.getParameter< std::string >("trackType")), + //muonsTag(parset.getParameter< edm::InputTag >("muonsTag")), + //vxtTag(parset.getParameter< edm::InputTag >("vxtTag")), + //useIPxy(parset.getUntrackedParameter< bool >("useIPxy", true)), + //useIPz(parset.getUntrackedParameter< bool >("useIPz", true)), + //selectionTags(parset.getParameter< std::vector >("selectionTags")), + //trackType(parset.getParameter< std::string >("trackType")), parset_(parset) { produces(); diff --git a/Validation/RecoMuon/python/associators_cff.py b/Validation/RecoMuon/python/associators_cff.py index 0e2d19471d3bf..c705f3d75cfab 100644 --- a/Validation/RecoMuon/python/associators_cff.py +++ b/Validation/RecoMuon/python/associators_cff.py @@ -32,6 +32,16 @@ extractedGlobalMuons.trackType = "globalTrack" extractedMuonTracks_seq = cms.Sequence( extractedGlobalMuons ) +#----------ME0Muon Collection Production for association by chi2 +me0muon = cms.EDProducer("ME0MuonTrackCollProducer", + me0MuonTag = cms.InputTag("me0SegmentMatching"), + selectionTags = cms.vstring('All'), + ) +#-------------------- +me0muonColl_seq = cms.Sequence( + me0muon + ) + # # Configuration for Seed track extractor #