Skip to content

Commit

Permalink
Merge pull request #4239 from cerati/mtv-add-eff-vs-dr
Browse files Browse the repository at this point in the history
MultiTrackValidator: add efficiency and fakerate vs deltaR
  • Loading branch information
ktf committed Jun 17, 2014
2 parents 6739b70 + ea1f241 commit 41e6c7d
Show file tree
Hide file tree
Showing 11 changed files with 145 additions and 17 deletions.
9 changes: 7 additions & 2 deletions Validation/RecoTrack/interface/MTVHistoProducerAlgo.h
Expand Up @@ -44,7 +44,7 @@ class MTVHistoProducerAlgo{
const TrackingParticle::Vector& momentumTP, const TrackingParticle::Point& vertexTP,
double dxy, double dz, int nSimHits,
const reco::Track* track,
int numVertices, double vertz)=0;
int numVertices, double vertz, double dR)=0;

virtual void fill_recoAssociated_simTrack_histos(int count,
const reco::GenParticle& tp,
Expand All @@ -63,7 +63,7 @@ class MTVHistoProducerAlgo{
int numVertices,
int tpbunchcrossing,
int nSimHits,
double sharedFraction)=0;
double sharedFraction, double dR)=0;

virtual void fill_dedx_recoTrack_histos(int count, edm::RefToBase<reco::Track>& trackref,const std::vector< edm::ValueMap<reco::DeDxData> >& v_dEdx)=0;
// virtual void fill_dedx_recoTrack_histos(reco::TrackRef trackref, std::vector< edm::ValueMap<reco::DeDxData> > v_dEdx)=0;
Expand Down Expand Up @@ -110,6 +110,11 @@ class MTVHistoProducerAlgo{
std::vector<int>& denominator,
std::string type);

void fillPlotFromPlots(MonitorElement* h,
TH1* numerator,
TH1* denominator,
std::string type);

void BinLogX(TH1*h);

private:
Expand Down
Expand Up @@ -45,7 +45,7 @@ class MTVHistoProducerAlgoForTracker: public MTVHistoProducerAlgo {
const TrackingParticle::Vector& momentumTP, const TrackingParticle::Point& vertexTP,
double dxy, double dz, int nSimHits,
const reco::Track* track,
int numVertices, double vertz);
int numVertices, double vertz, double dR);

void fill_recoAssociated_simTrack_histos(int count,
const reco::GenParticle& tp,
Expand All @@ -65,7 +65,8 @@ class MTVHistoProducerAlgoForTracker: public MTVHistoProducerAlgo {
int numVertices,
int tpbunchcrossing,
int nSimHits,
double sharedFraction);
double sharedFraction,
double dR);

void fill_dedx_recoTrack_histos(int count, edm::RefToBase<reco::Track>& trackref, const std::vector< edm::ValueMap<reco::DeDxData> >& v_dEdx);
// void fill_dedx_recoTrack_histos(reco::TrackRef trackref, std::vector< edm::ValueMap<reco::DeDxData> > v_dEdx);
Expand Down Expand Up @@ -139,6 +140,7 @@ class MTVHistoProducerAlgoForTracker: public MTVHistoProducerAlgo {
double minDz, maxDz; int nintDz;
double minVertpos, maxVertpos; int nintVertpos;
double minZpos, maxZpos; int nintZpos;
double mindr, maxdr; int nintdr;
double minDeDx, maxDeDx; int nintDeDx;
double minVertcount, maxVertcount; int nintVertcount;

Expand All @@ -163,6 +165,7 @@ class MTVHistoProducerAlgoForTracker: public MTVHistoProducerAlgo {
std::vector<MonitorElement*> h_effic_vs_dz, h_fake_vs_dz, h_recodz, h_assocdz, h_assoc2dz, h_simuldz, h_looperdz, h_misiddz, h_loopratedz, h_misidratedz;

std::vector<MonitorElement*> h_effic_vs_vertpos, h_effic_vs_zpos, h_assocvertpos, h_simulvertpos, h_assoczpos, h_simulzpos;
std::vector<MonitorElement*> h_effic_vs_dr, h_fakerate_vs_dr, h_assocdr, h_assoc2dr, h_simuldr, h_recodr;
std::vector<MonitorElement*> h_pt, h_eta, h_pullTheta,h_pullPhi,h_pullDxy,h_pullDz,h_pullQoverp;
std::vector<MonitorElement*> h_effic_vertcount_entire, h_fakerate_vertcount_entire, h_reco_vertcount_entire, h_assoc_vertcount_entire, h_assoc2_vertcount_entire, h_simul_vertcount_entire;
std::vector<MonitorElement*> h_effic_vertcount_barrel, h_fakerate_vertcount_barrel, h_reco_vertcount_barrel, h_assoc_vertcount_barrel, h_assoc2_vertcount_barrel, h_simul_vertcount_barrel;
Expand Down
3 changes: 3 additions & 0 deletions Validation/RecoTrack/interface/MultiTrackValidator.h
Expand Up @@ -46,8 +46,11 @@ class MultiTrackValidator : public thread_unsafe::DQMEDAnalyzer, protected Multi
//(i.e. "denominator" of the efficiency ratio)
TrackingParticleSelector tpSelector;
CosmicTrackingParticleSelector cosmictpSelector;
TrackingParticleSelector dRtpSelector;

edm::EDGetTokenT<SimHitTPAssociationProducer::SimHitTPAssociationList> _simHitTpMapTag;
edm::EDGetTokenT<edm::View<reco::Track> > labelTokenForDrCalculation;

};


Expand Down
59 changes: 56 additions & 3 deletions Validation/RecoTrack/plugins/MultiTrackValidator.cc
Expand Up @@ -31,6 +31,7 @@

#include "TMath.h"
#include <TF1.h>
#include "DataFormats/Math/interface/deltaR.h"

//#include <iostream>

Expand Down Expand Up @@ -75,11 +76,26 @@ MultiTrackValidator::MultiTrackValidator(const edm::ParameterSet& pset):MultiTra
pset.getParameter<bool>("chargedOnlyTP"),
pset.getParameter<std::vector<int> >("pdgIdTP"));


ParameterSet psetVsEta = psetForHistoProducerAlgo.getParameter<ParameterSet>("TpSelectorForEfficiencyVsEta");
dRtpSelector = TrackingParticleSelector(psetVsEta.getParameter<double>("ptMin"),
psetVsEta.getParameter<double>("minRapidity"),
psetVsEta.getParameter<double>("maxRapidity"),
psetVsEta.getParameter<double>("tip"),
psetVsEta.getParameter<double>("lip"),
psetVsEta.getParameter<int>("minHit"),
psetVsEta.getParameter<bool>("signalOnly"),
psetVsEta.getParameter<bool>("chargedOnly"),
psetVsEta.getParameter<bool>("stableOnly"),
psetVsEta.getParameter<std::vector<int> >("pdgId"));

useGsf = pset.getParameter<bool>("useGsf");
runStandalone = pset.getParameter<bool>("runStandalone");

_simHitTpMapTag = mayConsume<SimHitTPAssociationProducer::SimHitTPAssociationList>(pset.getParameter<edm::InputTag>("simHitTpMapTag"));

labelTokenForDrCalculation = consumes<edm::View<reco::Track> >(pset.getParameter<edm::InputTag>("trackCollectionForDrCalculation"));

if (!UseAssociators) {
associators.clear();
associators.push_back(assMapInput.label());
Expand Down Expand Up @@ -196,6 +212,28 @@ void MultiTrackValidator::analyze(const edm::Event& event, const edm::EventSetup
event.getByToken(label_tv,tvH);
TrackingVertexCollection tv = *tvH;

//calculate dR for TPs
std::vector<double> dR_tPCeff;
for (TrackingParticleCollection::size_type i=0; i<tPCeff.size(); i++){
TrackingParticleRef tpr(TPCollectionHeff, i);
TrackingParticle* tp=const_cast<TrackingParticle*>(tpr.get());
double dR = std::numeric_limits<double>::max();
if(dRtpSelector(*tp)) {//only for those needed for efficiency!
for (TrackingParticleCollection::size_type j=0; j<tPCeff.size(); j++){
if (i==j) continue;
TrackingParticleRef tpr2(TPCollectionHeff, j);
TrackingParticle* tp2=const_cast<TrackingParticle*>(tpr2.get());
if(! tpSelector(*tp2)) continue;//calculare dR wrt inclusive collection (also with PU, low pT, displaced)
double dR_tmp = reco::deltaR(*tp,*tp2);
if (dR_tmp<dR) dR=dR_tmp;
}
}
dR_tPCeff.push_back(dR);
}

edm::Handle<View<Track> > trackCollectionForDrCalculation;
event.getByToken(labelTokenForDrCalculation, trackCollectionForDrCalculation);

int w=0; //counter counting the number of sets of histograms
for (unsigned int ww=0;ww<associators.size();ww++){
for (unsigned int www=0;www<label.size();www++){
Expand Down Expand Up @@ -265,6 +303,7 @@ void MultiTrackValidator::analyze(const edm::Event& event, const edm::EventSetup
TrackingParticle::Point vertexTP;
double dxySim(0);
double dzSim(0);
double dR=dR_tPCeff[i];

//---------- THIS PART HAS TO BE CLEANED UP. THE PARAMETER DEFINER WAS NOT MEANT TO BE USED IN THIS WAY ----------
//If the TrackingParticle is collison like, get the momentum and vertex at production state
Expand Down Expand Up @@ -342,7 +381,7 @@ void MultiTrackValidator::analyze(const edm::Event& event, const edm::EventSetup
}
}

histoProducerAlgo_->fill_recoAssociated_simTrack_histos(w,*tp,momentumTP,vertexTP,dxySim,dzSim,nSimHits,matchedTrackPointer,puinfo.getPU_NumInteractions(), vtx_z_PU);
histoProducerAlgo_->fill_recoAssociated_simTrack_histos(w,*tp,momentumTP,vertexTP,dxySim,dzSim,nSimHits,matchedTrackPointer,puinfo.getPU_NumInteractions(), vtx_z_PU,dR);
sts++;
if (matchedTrackPointer) asts++;

Expand Down Expand Up @@ -390,6 +429,20 @@ void MultiTrackValidator::analyze(const edm::Event& event, const edm::EventSetup
}
//end dE/dx


//calculate dR for tracks
std::vector<double> dR_trk;
for(View<Track>::size_type i=0; i<trackCollection->size(); ++i){
RefToBase<Track> track(trackCollection, i);
double dR = std::numeric_limits<double>::max();
for(View<Track>::size_type j=0; j<trackCollectionForDrCalculation->size(); ++j){
RefToBase<Track> track2(trackCollectionForDrCalculation, j);
double dR_tmp = reco::deltaR(*track,*track2);
if (dR_tmp<dR && dR_tmp>std::numeric_limits<double>::min()) dR=dR_tmp;
}
dR_trk.push_back(dR);
}

for(View<Track>::size_type i=0; i<trackCollection->size(); ++i){

RefToBase<Track> track(trackCollection, i);
Expand Down Expand Up @@ -430,8 +483,8 @@ void MultiTrackValidator::analyze(const edm::Event& event, const edm::EventSetup
<< " NOT associated to any TrackingParticle" << "\n";
}


histoProducerAlgo_->fill_generic_recoTrack_histos(w,*track,bs.position(),isSimMatched,isSigSimMatched, isChargeMatched, numAssocRecoTracks, puinfo.getPU_NumInteractions(), tpbx, nSimHits, sharedFraction);
double dR=dR_trk[i];
histoProducerAlgo_->fill_generic_recoTrack_histos(w,*track,bs.position(),isSimMatched,isSigSimMatched, isChargeMatched, numAssocRecoTracks, puinfo.getPU_NumInteractions(), tpbx, nSimHits, sharedFraction,dR);

// dE/dx
// reco::TrackRef track2 = reco::TrackRef( trackCollection, i );
Expand Down
3 changes: 2 additions & 1 deletion Validation/RecoTrack/plugins/MultiTrackValidatorGenPs.cc
Expand Up @@ -345,7 +345,8 @@ void MultiTrackValidatorGenPs::analyze(const edm::Event& event, const edm::Event
}


histoProducerAlgo_->fill_generic_recoTrack_histos(w,*track,bs.position(),isGenMatched,isSigGenMatched, isChargeMatched, numAssocRecoTracks, puinfo.getPU_NumInteractions(), tpbx, nSimHits, sharedFraction);
double dR=0;//fxime
histoProducerAlgo_->fill_generic_recoTrack_histos(w,*track,bs.position(),isGenMatched,isSigGenMatched, isChargeMatched, numAssocRecoTracks, puinfo.getPU_NumInteractions(), tpbx, nSimHits, sharedFraction,dR);

// dE/dx
// reco::TrackRef track2 = reco::TrackRef( trackCollection, i );
Expand Down
6 changes: 4 additions & 2 deletions Validation/RecoTrack/plugins/TrackerSeedValidator.cc
Expand Up @@ -261,8 +261,9 @@ void TrackerSeedValidator::analyze(const edm::Event& event, const edm::EventSetu
matchedTrackPointer->setHitPattern(matchedSeedPointer->recHits().first,matchedSeedPointer->recHits().second);
}

double dR=0;//fixme
histoProducerAlgo_->fill_recoAssociated_simTrack_histos(w,*tp,tp->momentum(),tp->vertex(),dxySim,dzSim,nSimHits,
matchedTrackPointer,puinfo.getPU_NumInteractions(), vtx_z_PU);
matchedTrackPointer,puinfo.getPU_NumInteractions(), vtx_z_PU,dR);

sts++;
if (matchedTrackPointer) asts++;
Expand Down Expand Up @@ -345,9 +346,10 @@ void TrackerSeedValidator::analyze(const edm::Event& event, const edm::EventSetu
edm::LogVerbatim("SeedValidator") << "TrajectorySeed #" << rT << " NOT associated to any TrackingParticle" << "\n";
}

double dR = 0.;//fixme
histoProducerAlgo_->fill_generic_recoTrack_histos(w,*trackFromSeed,bs.position(),isSimMatched,isSigSimMatched,
isChargeMatched, numAssocSeeds, puinfo.getPU_NumInteractions(),
tpbx, nSimHits, sharedFraction);
tpbx, nSimHits, sharedFraction, dR);

//Fill other histos
try{
Expand Down
Expand Up @@ -69,6 +69,11 @@
minZpos = cms.double(-30),
maxZpos = cms.double(30),
nintZpos = cms.int32(60),
#
# dR
mindr = cms.double(0.001),
maxdr = cms.double(1),
nintdr = cms.int32(100),

# Pileup vertices
minVertcount = cms.double(-0.5),
Expand Down
3 changes: 3 additions & 0 deletions Validation/RecoTrack/python/MultiTrackValidator_cfi.py
Expand Up @@ -57,4 +57,7 @@
### output configuration
dirName = cms.string('Tracking/Track/'),
outputFile = cms.string(''),

### for fake rate vs dR ###
trackCollectionForDrCalculation = cms.InputTag("generalTracks")
)
2 changes: 2 additions & 0 deletions Validation/RecoTrack/python/PostProcessorTracker_cfi.py
Expand Up @@ -23,6 +23,7 @@
"chargeMisIdRate_dz 'Charge MisID Rate vs Dz' num_chargemisid_versus_dz num_reco_dz",
"effic_vs_vertpos 'Efficiency vs vertpos' num_assoc(simToReco)_vertpos num_simul_vertpos",
"effic_vs_zpos 'Efficiency vs zpos' num_assoc(simToReco)_zpos num_simul_zpos",
"effic_vs_dr 'Efficiency vs dr' num_assoc(simToReco)_dr num_simul_dr",
"effic_vertcount_barrel 'efficiency in barrel vs N of pileup vertices' num_assoc(simToReco)_vertcount_barrel num_simul_vertcount_barrel",
"effic_vertcount_fwdpos 'efficiency in endcap(+) vs N of pileup vertices' num_assoc(simToReco)_vertcount_fwdpos num_simul_vertcount_fwdpos",
"effic_vertcount_fwdneg 'efficiency in endcap(-) vs N of pileup vertices' num_assoc(simToReco)_vertcount_fwdneg num_simul_vertcount_fwdneg",
Expand All @@ -35,6 +36,7 @@
"fakerate_vs_phi 'Fake rate vs phi' num_assoc(recoToSim)_phi num_reco_phi fake",
"fakerate_vs_dxy 'Fake rate vs dxy' num_assoc(recoToSim)_dxy num_reco_dxy fake",
"fakerate_vs_dz 'Fake rate vs dz' num_assoc(recoToSim)_dz num_reco_dz fake",
"fakerate_vs_dr 'Fake rate vs dr' num_assoc(recoToSim)_dr num_reco_dr fake",
"fakerate_vertcount_barrel 'fake rate in barrel vs N of pileup vertices' num_assoc(recoToSim)_vertcount_barrel num_reco_vertcount_barrel fake",
"fakerate_vertcount_fwdpos 'fake rate in endcap(+) vs N of pileup vertices' num_assoc(recoToSim)_vertcount_fwdpos num_reco_vertcount_fwdpos fake",
"fakerate_vertcount_fwdneg 'fake rate in endcap(-) vs N of pileup vertices' num_assoc(recoToSim)_vertcount_fwdneg num_reco_vertcount_fwdneg fake"
Expand Down
26 changes: 24 additions & 2 deletions Validation/RecoTrack/src/MTVHistoProducerAlgo.cc
Expand Up @@ -44,8 +44,30 @@ void MTVHistoProducerAlgo::fillPlotFromVectors(MonitorElement* h,
}
}



void MTVHistoProducerAlgo::fillPlotFromPlots(MonitorElement* h,
TH1* numerator,
TH1* denominator,
std::string type){
assert(h->getNbinsX()==numerator->GetNbinsX());
assert(h->getNbinsX()==denominator->GetNbinsX());
double value,err;
for (int bin=1;bin<=h->getNbinsX();++bin) {
if (denominator->GetBinContent(bin)!=0) {
if (type=="effic"){
value = ((double) numerator->GetBinContent(bin))/((double) denominator->GetBinContent(bin));
err = sqrt( value*(1.-value)/((double) denominator->GetBinContent(bin)) );
} else if (type=="fakerate"){
value = 1.-((double) numerator->GetBinContent(bin))/((double) denominator->GetBinContent(bin));
err = sqrt( value*(1.-value)/(double)denominator->GetBinContent(bin) );
} else if (type=="pileup"){
value = ((double) numerator->GetBinContent(bin))/((double) denominator->GetBinContent(bin));
err = sqrt( value*(1.+value)/(double) denominator->GetBinContent(bin) );
} else return;
h->setBinContent(bin,value);
h->setBinError(bin,err);
}
}
}

void MTVHistoProducerAlgo::BinLogX(TH1*h){
TAxis *axis = h->GetXaxis();
Expand Down

0 comments on commit 41e6c7d

Please sign in to comment.