CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/DQM/HLTEvF/plugins/HLTMonBTagIPSource.cc

Go to the documentation of this file.
00001 
00011 #include <vector>
00012 #include <string>
00013 #include <algorithm>
00014 
00015 #include "DQMServices/Core/interface/DQMStore.h"
00016 #include "DQMServices/Core/interface/MonitorElement.h"
00017 #include "FWCore/Framework/interface/Event.h"
00018 #include "FWCore/Framework/interface/EventSetup.h"
00019 #include "FWCore/Framework/interface/LuminosityBlock.h"
00020 #include "FWCore/Framework/interface/Run.h"
00021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00022 #include "FWCore/Utilities/interface/Exception.h"
00023 #include "DataFormats/Common/interface/Handle.h"
00024 #include "DataFormats/Common/interface/View.h"
00025 #include "DataFormats/Common/interface/TriggerResults.h"
00026 #include "DataFormats/BTauReco/interface/TrackIPTagInfo.h"
00027 #include "DataFormats/BTauReco/interface/JetTag.h"
00028 #include "DataFormats/TrackReco/interface/Track.h"
00029 #include "DataFormats/JetReco/interface/Jet.h"
00030 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
00031 #include "HLTMonBTagIPSource.h"
00032 
00033 HLTMonBTagIPSource::HLTMonBTagIPSource(const edm::ParameterSet & config) :
00034   m_L1Filter(       config.getParameter<edm::InputTag>("L1Filter") ),
00035   m_L2Filter(       config.getParameter<edm::InputTag>("L2Filter") ),
00036   m_L25Filter(      config.getParameter<edm::InputTag>("L25Filter") ),
00037   m_L3Filter(       config.getParameter<edm::InputTag>("L3Filter") ),
00038   m_L2Jets(         config.getParameter<edm::InputTag>("L2Jets") ),
00039   m_L25TagInfo(     config.getParameter<edm::InputTag>("L25TagInfo") ),
00040   m_L25JetTags(     config.getParameter<edm::InputTag>("L25JetTags") ),
00041   m_L3TagInfo(      config.getParameter<edm::InputTag>("L3TagInfo") ),
00042   m_L3JetTags(      config.getParameter<edm::InputTag>("L3JetTags") ),
00043   m_triggerResults( config.getParameter<edm::InputTag>("triggerResults") ),
00044   m_processName(    config.getParameter<std::string>("processName") ),
00045   m_pathName(       config.getParameter<std::string>("pathName") ),
00046   m_monitorName(    config.getParameter<std::string>("monitorName" ) ),
00047   m_outputFile(     config.getUntrackedParameter<std::string>("outputFile", "HLTBJetDQM.root") ),
00048   m_storeROOT(      config.getUntrackedParameter<bool>("storeROOT", false) ),
00049   m_size(           config.getParameter<unsigned int>("interestingJets") ),
00050   m_dbe(),
00051   m_init(           false ),
00052   m_pathIndex(      (unsigned int) -1 ),
00053   m_L1FilterIndex(  (unsigned int) -1 ),
00054   m_L2FilterIndex(  (unsigned int) -1 ),
00055   m_L25FilterIndex( (unsigned int) -1 ),
00056   m_L3FilterIndex(  (unsigned int) -1 ),
00057 
00058   // MonitorElement's (plots) filled by the source
00059   m_plotRates(0),
00060 
00061   m_plotL2JetsEnergy(0),
00062   m_plotL2JetsET(0),
00063   m_plotL2JetsEta(0),
00064   m_plotL2JetsPhi(0),
00065   m_plotL2JetsEtaPhi(0),
00066   m_plotL2JetsEtaET(0),
00067   m_plotL25JetsEnergy(0),
00068   m_plotL25JetsET(0),
00069   m_plotL25JetsEta(0),
00070   m_plotL25JetsPhi(0),
00071   m_plotL25JetsEtaPhi(0),
00072   m_plotL25JetsEtaET(0),
00073   m_plotL25TrackMultiplicity(0),
00074   m_plotL25TrackHits(0),
00075   m_plotL25TrackChi2(0),
00076   m_plotL25TrackEtaPhi(0),
00077   m_plotL25TrackEtaPT(0),
00078   m_plotL25IP2ndTrack2d(0),
00079   m_plotL25IP2ndTrack2dSig(0),
00080   m_plotL25IP2ndTrack3d(0),
00081   m_plotL25IP2ndTrack3dSig(0),
00082   m_plotL25IP3ndTrack2d(0),
00083   m_plotL25IP3ndTrack2dSig(0),
00084   m_plotL25IP3ndTrack3d(0),
00085   m_plotL25IP3ndTrack3dSig(0),
00086   m_plotL25Discriminator(0),
00087   m_plotL3JetsEnergy(0),
00088   m_plotL3JetsET(0),
00089   m_plotL3JetsEta(0),
00090   m_plotL3JetsPhi(0),
00091   m_plotL3JetsEtaPhi(0),
00092   m_plotL3JetsEtaET(0),
00093   m_plotL3TrackMultiplicity(0),
00094   m_plotL3TrackHits(0),
00095   m_plotL3TrackChi2(0),
00096   m_plotL3TrackEtaPhi(0),
00097   m_plotL3TrackEtaPT(0),
00098   m_plotL3IP2ndTrack2d(0),
00099   m_plotL3IP2ndTrack2dSig(0),
00100   m_plotL3IP2ndTrack3d(0),
00101   m_plotL3IP2ndTrack3dSig(0),
00102   m_plotL3IP3ndTrack2d(0),
00103   m_plotL3IP3ndTrack2dSig(0),
00104   m_plotL3IP3ndTrack3d(0),
00105   m_plotL3IP3ndTrack3dSig(0),
00106   m_plotL3Discriminator(0)
00107 {
00108 }
00109 
00110 HLTMonBTagIPSource::~HLTMonBTagIPSource(void) {
00111 }
00112 
00113 void HLTMonBTagIPSource::beginJob() {
00114   if (not m_dbe.isAvailable())
00115     return;
00116 
00117   m_dbe->setVerbose(0);
00118   m_dbe->setCurrentFolder(m_monitorName + "/" + m_pathName);
00119 
00120   m_plotRates                       = book("Rates",                  "Rates",                              6,  0.,     6);
00121 
00122   m_plotL2JetsEnergy                = book("L2_jet_energy",          "L2 jet energy",                    300,   0.,  300.,  "GeV");
00123   m_plotL2JetsET                    = book("L2_jet_eT",              "L2 jet eT",                        300,   0.,  300.,  "GeV");
00124   m_plotL2JetsEta                   = book("L2_jet_eta",             "L2 jet eta",                        60,  -3.0,   3.0, "#eta");
00125   m_plotL2JetsPhi                   = book("L2_jet_phi",             "L2 jet phi",                        64,  -3.2,   3.2, "#phi");
00126   m_plotL2JetsEtaPhi                = book("L2_jet_eta_phi",         "L2 jet eta vs. phi",                60,  -3.0,   3.0,  64, -3.2,   3.2, "#eta", "#phi");
00127   m_plotL2JetsEtaET                 = book("L2_jet_eta_et",          "L2 jet eta vs. eT",                 60,  -3.0,   3.0, 300,  0.,  300.,  "#eta", "GeV");
00128   m_plotL25JetsEnergy               = book("L25_jet_energy",         "L2.5 jet Energy",                  300,   0.,  300.,  "GeV");
00129   m_plotL25JetsET                   = book("L25_jet_eT",             "L2.5 jet ET",                      300,   0.,  300.,  "GeV");
00130   m_plotL25JetsEta                  = book("L25_jet_eta",            "L2.5 jet eta",                      60,  -3.0,   3.0, "#eta");
00131   m_plotL25JetsPhi                  = book("L25_jet_phi",            "L2.5 jet phi",                      64,  -3.2,   3.2, "#phi");
00132   m_plotL25JetsEtaPhi               = book("L25_jet_eta_phi",        "L2.5 jet eta vs. phi",              60,  -3.0,   3.0,  64, -3.2,   3.2, "#eta", "#phi");
00133   m_plotL25JetsEtaET                = book("L25_jet_eta_et",         "L2.5 jet eta vs. eT",               60,  -3.0,   3.0, 300,  0.,  300.,  "#eta", "GeV");
00134   m_plotL25TrackMultiplicity        = book("L25_track_multiplicity", "L2.5 pixel tracks multiplicity",    25,   0,    25);
00135   m_plotL25TrackHits                = book("L25_track_hits",         "L2.5 pixel tracks n. of hits",       5,   0,     5);
00136   m_plotL25TrackChi2                = book("L25_track_chi2",         "L2.5 pixel tracks chi2/DoF",        20,   0.,   20.,  "#chi^2/DoF");
00137   m_plotL25TrackEtaPhi              = book("L25_track_eta_phi",      "L2.5 pixel tracks eta vs. phi",     60,  -3.0,   3.0,  64, -3.2,   3.2, "#eta", "#phi");
00138   m_plotL25TrackEtaPT               = book("L25_track_eta_pt",       "L2.5 pixel tracks eta vs. pT",      60,  -3.0,   3.0,  50,  0.,   50.,  "#eta", "GeV");
00139   m_plotL25IP2ndTrack2d             = book("L25_IP_2ndTrack_2d",     "L2.5 2nd pixel track 2D IP",        25,  -0.05, 0.20, "cm");
00140   m_plotL25IP2ndTrack2dSig          = book("L25_IP_2ndTrack_2dSig",  "L2.5 2nd pixel track 2D SIP",       80, -30.,   50.);
00141   m_plotL25IP2ndTrack3d             = book("L25_IP_2ndTrack_3d",     "L2.5 2nd pixel track 3D IP",        60,  -0.20, 1.00, "cm");
00142   m_plotL25IP2ndTrack3dSig          = book("L25_IP_2ndTrack_3dSig",  "L2.5 2nd pixel track 3D SIP",       80, -30.,   50.);
00143   m_plotL25IP3ndTrack2d             = book("L25_IP_3ndTrack_2d",     "L2.5 3rd pixel track 2D IP",        25,  -0.05, 0.20, "cm");
00144   m_plotL25IP3ndTrack2dSig          = book("L25_IP_3ndTrack_2dSig",  "L2.5 3rd pixel track 2D SIP",       80, -30.,   50.);
00145   m_plotL25IP3ndTrack3d             = book("L25_IP_3ndTrack_3d",     "L2.5 3rd pixel track 3D IP",        60,  -0.20, 1.00, "cm");
00146   m_plotL25IP3ndTrack3dSig          = book("L25_IP_3ndTrack_3dSig",  "L2.5 3rd pixel track 3D SIP",       80, -30.,   50.);
00147   m_plotL25Discriminator            = book("L25_discriminator",      "L2.5 b-tag discriminator",          80, -30.,   50.);
00148   m_plotL3JetsEnergy                = book("L3_jet_energy",          "L3 jet Energy",                    300,   0.,  300.,  "GeV");
00149   m_plotL3JetsET                    = book("L3_jet_eT",              "L3 jet ET",                        300,   0.,  300.,  "GeV");
00150   m_plotL3JetsEta                   = book("L3_jet_eta",             "L3 jet eta",                        60,  -3.0,   3.0, "#eta");
00151   m_plotL3JetsPhi                   = book("L3_jet_phi",             "L3 jet phi",                        64,  -3.2,   3.2, "#phi");
00152   m_plotL3JetsEtaPhi                = book("L3_jet_eta_phi",         "L3 jet eta vs. phi",                60,  -3.0,   3.0,  64, -3.2,   3.2, "#eta", "#phi");
00153   m_plotL3JetsEtaET                 = book("L3_jet_eta_et",          "L3 jet eta vs. eT",                 60,  -3.0,   3.0, 300,  0.,  300.,  "#eta", "GeV");
00154   m_plotL3TrackMultiplicity         = book("L3_track_multiplicity",  "L3 tracks multiplicity",            25,   0,    25);
00155   m_plotL3TrackHits                 = book("L3_track_hits",          "L3 tracks n. of hits",              20,   0,    20);
00156   m_plotL3TrackChi2                 = book("L3_track_chi2",          "L3 tracks chi2/DoF",                20,   0.,   20.,  "#chi^2/DoF");
00157   m_plotL3TrackEtaPhi               = book("L3_track_eta_phi",       "L3 tracks eta vs. phi",             60,  -3.0,   3.0,  64, -3.2,   3.2, "#eta", "#phi");
00158   m_plotL3TrackEtaPT                = book("L3_track_eta_pt",        "L3 tracks eta vs. pT",              60,  -3.0,   3.0,  50,  0.,   50.,  "#eta", "GeV");
00159   m_plotL3IP2ndTrack2d              = book("L3_IP_2ndTrack_2d",      "L3 2nd track 2D IP",                25,  -0.05, 0.20, "cm");
00160   m_plotL3IP2ndTrack2dSig           = book("L3_IP_2ndTrack_2dSig",   "L3 2nd track 2D SIP",               80, -30.,   50.);
00161   m_plotL3IP2ndTrack3d              = book("L3_IP_2ndTrack_3d",      "L3 2nd track 3D IP",                60,  -0.20, 1.00, "cm");
00162   m_plotL3IP2ndTrack3dSig           = book("L3_IP_2ndTrack_3dSig",   "L3 2nd track 3D SIP",               80, -30.,   50.);
00163   m_plotL3IP3ndTrack2d              = book("L3_IP_3ndTrack_2d",      "L3 3rd track 2D IP",                25,  -0.05, 0.20, "cm");
00164   m_plotL3IP3ndTrack2dSig           = book("L3_IP_3ndTrack_2dSig",   "L3 3rd track 2D SIP",               80, -30.,   50.);
00165   m_plotL3IP3ndTrack3d              = book("L3_IP_3ndTrack_3d",      "L3 3rd track 3D IP",                60,  -0.20, 1.00, "cm");
00166   m_plotL3IP3ndTrack3dSig           = book("L3_IP_3ndTrack_3dSig",   "L3 3rd track 3D SIP",               80, -30.,   50.);
00167   m_plotL3Discriminator             = book("L3_discriminator",       "L3 b-tag discriminator",            80, -30.,   50.);
00168 }
00169 
00170 void HLTMonBTagIPSource::endJob() { 
00171   if (m_dbe.isAvailable() and m_storeROOT)
00172     m_dbe->save(m_outputFile);
00173 }
00174 
00175 void HLTMonBTagIPSource::beginRun(const edm::Run & run, const edm::EventSetup & setup) {
00176 
00177   HLTConfigProvider configProvider;
00178 
00179   bool changed = false;
00180   if (not configProvider.init(run, setup, m_processName, changed))
00181   {
00182     edm::LogWarning("ConfigurationError") << "process name \"" << m_processName << "\" is not valid.";
00183     m_init = false;
00184     return;
00185   }
00186 
00187   m_pathIndex = configProvider.triggerIndex( m_pathName );
00188   if (m_pathIndex == configProvider.size())
00189   {
00190     edm::LogWarning("ConfigurationError") << "trigger name \"" << m_pathName << "\" is not valid.";
00191     m_init = false;
00192     return;
00193   }
00194 
00195   m_init = true;
00196 
00197   // if their call fails, these will be set to one after the last valid module for their path
00198   // so they will never be "passed"
00199   unsigned int size = configProvider.size( m_pathIndex );
00200 
00201   m_L1FilterIndex  = configProvider.moduleIndex( m_pathIndex, m_L1Filter.encode()  );
00202   if (m_L1FilterIndex == size)
00203     edm::LogWarning("ConfigurationError") << "L1 filter \"" << m_L1Filter << "\" is not valid.";
00204 
00205   m_L2FilterIndex  = configProvider.moduleIndex( m_pathIndex, m_L2Filter.encode()  );
00206   if (m_L2FilterIndex == size)
00207     edm::LogWarning("ConfigurationError") << "L2 filter \"" << m_L2Filter << "\" is not valid.";
00208 
00209   m_L25FilterIndex = configProvider.moduleIndex( m_pathIndex, m_L25Filter.encode() );
00210   if (m_L25FilterIndex == size)
00211     edm::LogWarning("ConfigurationError") << "L2.5 filter \"" << m_L25Filter << "\" is not valid.";
00212 
00213   m_L3FilterIndex  = configProvider.moduleIndex( m_pathIndex, m_L3Filter.encode()  );
00214   if (m_L3FilterIndex == size)
00215     edm::LogWarning("ConfigurationError") << "L3 filter \"" << m_L3Filter << "\" is not valid.";
00216 }
00217 
00218 void HLTMonBTagIPSource::endRun(const edm::Run & run, const edm::EventSetup & setup) {
00219 }
00220 
00221 void HLTMonBTagIPSource::beginLuminosityBlock(const edm::LuminosityBlock & lumi, const edm::EventSetup & setup) {
00222 }
00223 
00224 void HLTMonBTagIPSource::endLuminosityBlock(const edm::LuminosityBlock & lumi, const edm::EventSetup & setup) {
00225 }
00226 
00227 void HLTMonBTagIPSource::analyze(const edm::Event & event, const edm::EventSetup & setup) {
00228   if (not m_dbe.isAvailable())
00229     return;
00230 
00231   if (not m_init)
00232     return;
00233 
00234   edm::Handle<edm::TriggerResults> h_triggerResults;
00235   edm::Handle<edm::View<reco::Jet> >          h_L2Jets;
00236   edm::Handle<reco::TrackIPTagInfoCollection> h_L25TagInfo;
00237   edm::Handle<reco::JetTagCollection>         h_L25JetTags;
00238   edm::Handle<reco::TrackIPTagInfoCollection> h_L3TagInfo;
00239   edm::Handle<reco::JetTagCollection>         h_L3JetTags;
00240   
00241   event.getByLabel(m_triggerResults, h_triggerResults);
00242   event.getByLabel(m_L2Jets,     h_L2Jets);
00243   event.getByLabel(m_L25TagInfo, h_L25TagInfo);
00244   event.getByLabel(m_L25JetTags, h_L25JetTags);
00245   event.getByLabel(m_L3TagInfo,  h_L3TagInfo);
00246   event.getByLabel(m_L3JetTags,  h_L3JetTags);
00247 
00248   // check if this path passed the L1, L2, L2.5 and L3 filters
00249   bool         wasrun = false;
00250   unsigned int latest = 0;
00251   bool         accept = false;
00252   if (h_triggerResults.isValid()) {
00253     wasrun = h_triggerResults->wasrun( m_pathIndex );
00254     latest = h_triggerResults->index(  m_pathIndex );
00255     accept = h_triggerResults->accept( m_pathIndex );
00256   }
00257   if (wasrun)
00258     m_plotRates->Fill( 0. );    // path was run
00259   if (latest > m_L1FilterIndex)
00260     m_plotRates->Fill( 1. );    // L1 accepted
00261   if (latest > m_L2FilterIndex)
00262     m_plotRates->Fill( 2. );    // L2 accepted  
00263   if (latest > m_L25FilterIndex)
00264     m_plotRates->Fill( 3. );    // L2.5 accepted  
00265   if (latest > m_L3FilterIndex)
00266     m_plotRates->Fill( 4. );    // L3 accepted  
00267   if (accept)
00268     m_plotRates->Fill( 5. );    // HLT accepted
00269 
00270   if ((latest > m_L1FilterIndex) and h_L2Jets.isValid()) {
00271     unsigned int size = std::min((unsigned int) h_L2Jets->size(), m_size);
00272     for (unsigned int i = 0; i < size; ++i) {
00273       const reco::Jet & jet = (*h_L2Jets)[i];
00274       m_plotL2JetsEnergy->Fill( jet.energy() );
00275       m_plotL2JetsET->Fill(     jet.et() );
00276       m_plotL2JetsEta->Fill(    jet.eta() );
00277       m_plotL2JetsPhi->Fill(    jet.phi() );
00278       m_plotL2JetsEtaPhi->Fill( jet.eta(), jet.phi() );
00279       m_plotL2JetsEtaET->Fill(  jet.eta(), jet.et() );
00280     }
00281   }
00282   if ((latest > m_L2FilterIndex) and h_L25TagInfo.isValid() and h_L25JetTags.isValid()) {
00283     unsigned int size = std::min((unsigned int) h_L25TagInfo->size(), m_size);
00284     for (unsigned int i = 0; i < size; ++i) {
00285       const reco::TrackIPTagInfo & info   = (*h_L25TagInfo)[i];
00286       const reco::Jet & jet = * info.jet();
00287       const reco::TrackRefVector & tracks = info.selectedTracks();
00288       const std::vector<reco::TrackIPTagInfo::TrackIPData> & data = info.impactParameterData();
00289       const reco::JetTag & tag = (*h_L25JetTags)[info.jet().key()];
00290       m_plotL25JetsEnergy->Fill( jet.energy() );
00291       m_plotL25JetsET->Fill(     jet.et() );
00292       m_plotL25JetsEta->Fill(    jet.eta() );
00293       m_plotL25JetsPhi->Fill(    jet.phi() );
00294       m_plotL25JetsEtaPhi->Fill( jet.eta(), jet.phi() );
00295       m_plotL25JetsEtaET->Fill(  jet.eta(), jet.et() );
00296       m_plotL25TrackMultiplicity->Fill( tracks.size() );
00297       for (unsigned int t = 0; t < tracks.size(); ++t) {
00298         m_plotL25TrackHits->Fill(   tracks[t]->numberOfValidHits() );
00299         m_plotL25TrackChi2->Fill(   tracks[t]->normalizedChi2() );
00300         m_plotL25TrackEtaPhi->Fill( tracks[t]->eta(), tracks[t]->phi() );
00301         m_plotL25TrackEtaPT->Fill(  tracks[t]->eta(), tracks[t]->pt() );
00302       }
00303       std::vector<size_t> indicesBy2d = info.sortedIndexes(reco::TrackIPTagInfo::IP2DSig);
00304       if (indicesBy2d.size() >= 2) {
00305         m_plotL25IP2ndTrack2d->Fill(    data[indicesBy2d[1]].ip2d.value() );
00306         m_plotL25IP2ndTrack2dSig->Fill( data[indicesBy2d[1]].ip2d.significance() );
00307       }
00308       if (indicesBy2d.size() >= 3) {
00309         m_plotL25IP3ndTrack2d->Fill(    data[indicesBy2d[2]].ip2d.value() );
00310         m_plotL25IP3ndTrack2dSig->Fill( data[indicesBy2d[2]].ip2d.significance() );
00311       }
00312       std::vector<size_t> indicesBy3d = info.sortedIndexes(reco::TrackIPTagInfo::IP3DSig);
00313       if (indicesBy3d.size() >= 2) {
00314         m_plotL25IP2ndTrack3d->Fill(    data[indicesBy3d[1]].ip3d.value() );
00315         m_plotL25IP2ndTrack3dSig->Fill( data[indicesBy3d[1]].ip3d.significance() );
00316       }
00317       if (indicesBy3d.size() >= 3) {
00318         m_plotL25IP3ndTrack3d->Fill(    data[indicesBy3d[2]].ip3d.value() );
00319         m_plotL25IP3ndTrack3dSig->Fill( data[indicesBy3d[2]].ip3d.significance() );
00320       }
00321       m_plotL25Discriminator->Fill( tag.second );
00322     }
00323   }
00324   if ((latest > m_L25FilterIndex) and h_L3TagInfo.isValid() and h_L3JetTags.isValid()) {
00325     unsigned int size = std::min((unsigned int) h_L3TagInfo->size(), m_size);
00326     for (unsigned int i = 0; i < size; ++i) {
00327       const reco::TrackIPTagInfo & info   = (*h_L3TagInfo)[i];
00328       const reco::Jet & jet = * info.jet();
00329       const reco::TrackRefVector & tracks = info.selectedTracks();
00330       const std::vector<reco::TrackIPTagInfo::TrackIPData> & data = info.impactParameterData();
00331       const reco::JetTag & tag = (*h_L3JetTags)[info.jet().key()];
00332       m_plotL3JetsEnergy->Fill( jet.energy() );
00333       m_plotL3JetsET->Fill(     jet.et() );
00334       m_plotL3JetsEta->Fill(    jet.eta() );
00335       m_plotL3JetsPhi->Fill(    jet.phi() );
00336       m_plotL3JetsEtaPhi->Fill( jet.eta(), jet.phi() );
00337       m_plotL3JetsEtaET->Fill(  jet.eta(), jet.et() );
00338       m_plotL3TrackMultiplicity->Fill( tracks.size() );
00339       for (unsigned int t = 0; t < tracks.size(); ++t) {
00340         m_plotL3TrackHits->Fill(   tracks[t]->numberOfValidHits() );
00341         m_plotL3TrackChi2->Fill(   tracks[t]->normalizedChi2() );
00342         m_plotL3TrackEtaPhi->Fill( tracks[t]->eta(), tracks[t]->phi() );
00343         m_plotL3TrackEtaPT->Fill(  tracks[t]->eta(), tracks[t]->pt() );
00344       }
00345       std::vector<size_t> indicesBy2d = info.sortedIndexes(reco::TrackIPTagInfo::IP2DSig);
00346       if (indicesBy2d.size() >= 2) {
00347         m_plotL3IP2ndTrack2d->Fill(    data[indicesBy2d[1]].ip2d.value() );
00348         m_plotL3IP2ndTrack2dSig->Fill( data[indicesBy2d[1]].ip2d.significance() );
00349       }
00350       if (indicesBy2d.size() >= 3) {
00351         m_plotL3IP3ndTrack2d->Fill(    data[indicesBy2d[2]].ip2d.value() );
00352         m_plotL3IP3ndTrack2dSig->Fill( data[indicesBy2d[2]].ip2d.significance() );
00353       }
00354       std::vector<size_t> indicesBy3d = info.sortedIndexes(reco::TrackIPTagInfo::IP3DSig);
00355       if (indicesBy3d.size() >= 2) {
00356         m_plotL3IP2ndTrack3d->Fill(    data[indicesBy3d[1]].ip3d.value() );
00357         m_plotL3IP2ndTrack3dSig->Fill( data[indicesBy3d[1]].ip3d.significance() );
00358       }
00359       if (indicesBy3d.size() >= 3) {
00360         m_plotL3IP3ndTrack3d->Fill(    data[indicesBy3d[2]].ip3d.value() );
00361         m_plotL3IP3ndTrack3dSig->Fill( data[indicesBy3d[2]].ip3d.significance() );
00362       }
00363       m_plotL3Discriminator->Fill( tag.second );
00364     }
00365   }
00366 }
00367 
00368 MonitorElement * HLTMonBTagIPSource::book(const std::string & name, const std::string & title , int x_bins, double x_min, double x_max, const char * x_axis) {
00369   MonitorElement * element = m_dbe->book1D(name, title, x_bins, x_min, x_max);
00370   if (x_axis)
00371     element->setAxisTitle(x_axis, 1);
00372   return element;
00373 }
00374 
00375 MonitorElement * HLTMonBTagIPSource::book(const std::string & name, const std::string & title , int x_bins, double x_min, double x_max, int y_bins, double y_min, double y_max, const char * x_axis, const char * y_axis) {
00376   MonitorElement * element = m_dbe->book2D(name, title, x_bins, x_min, x_max, y_bins, y_min, y_max);
00377   if (x_axis)
00378     element->setAxisTitle(x_axis, 1);
00379   if (y_axis)
00380     element->setAxisTitle(y_axis, 2);
00381   return element;
00382 }
00383 
00384 // register as a framework plugin
00385 #include "FWCore/Framework/interface/MakerMacros.h"
00386 DEFINE_FWK_MODULE(HLTMonBTagIPSource);