CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_2_9_HLT1_bphpatch4/src/DQM/Physics/src/TopSingleLeptonDQM.cc

Go to the documentation of this file.
00001 #include "JetMETCorrections/Objects/interface/JetCorrectionsRecord.h"
00002 #include "DataFormats/JetReco/interface/CaloJet.h"
00003 #include "DataFormats/BTauReco/interface/JetTag.h"
00004 #include "DataFormats/JetReco/interface/PFJet.h"
00005 #include "DQM/Physics/src/TopSingleLeptonDQM.h"
00006 #include "DataFormats/Math/interface/deltaR.h"
00007 
00008 namespace TopSingleLepton {
00009 
00010   // maximal number of leading jets 
00011   // to be used for top mass estimate
00012   static const unsigned int MAXJETS = 4;
00013   // nominal mass of the W boson to 
00014   // be used for the top mass estimate
00015   static const double WMASS = 80.4;
00016 
00017   MonitorEnsemble::MonitorEnsemble(const char* label, const edm::ParameterSet& cfg) : 
00018     label_(label), elecIso_(0), elecSelect_(0), muonIso_(0), muonSelect_(0), jetIDSelect_(0), includeBTag_(false), lowerEdge_(-1.), upperEdge_(-1.), logged_(0)
00019   {
00020     // sources have to be given; this PSet is not optional
00021     edm::ParameterSet sources=cfg.getParameter<edm::ParameterSet>("sources");
00022     muons_= sources.getParameter<edm::InputTag>("muons");
00023     elecs_= sources.getParameter<edm::InputTag>("elecs");
00024     jets_ = sources.getParameter<edm::InputTag>("jets" );
00025     mets_ = sources.getParameter<std::vector<edm::InputTag> >("mets" );
00026 
00027     // electronExtras are optional; they may be omitted or 
00028     // empty
00029     if( cfg.existsAs<edm::ParameterSet>("elecExtras") ){
00030       edm::ParameterSet elecExtras=cfg.getParameter<edm::ParameterSet>("elecExtras");
00031       // select is optional; in case it's not found no
00032       // selection will be applied
00033       if( elecExtras.existsAs<std::string>("select") ){
00034         elecSelect_= new StringCutObjectSelector<reco::GsfElectron>(elecExtras.getParameter<std::string>("select"));
00035       }
00036       // isolation is optional; in case it's not found no
00037       // isolation will be applied
00038       if( elecExtras.existsAs<std::string>("isolation") ){
00039         elecIso_= new StringCutObjectSelector<reco::GsfElectron>(elecExtras.getParameter<std::string>("isolation"));
00040       }
00041       // electronId is optional; in case it's not found the 
00042       // InputTag will remain empty
00043       if( elecExtras.existsAs<edm::ParameterSet>("electronId") ){
00044         edm::ParameterSet elecId=elecExtras.getParameter<edm::ParameterSet>("electronId");
00045         electronId_= elecId.getParameter<edm::InputTag>("src");
00046         eidPattern_= elecId.getParameter<int>("pattern");
00047       }
00048     }
00049 
00050     // muonExtras are optional; they may be omitted or empty
00051     if( cfg.existsAs<edm::ParameterSet>("muonExtras") ){
00052       edm::ParameterSet muonExtras=cfg.getParameter<edm::ParameterSet>("muonExtras");
00053       // select is optional; in case it's not found no
00054       // selection will be applied
00055       if( muonExtras.existsAs<std::string>("select") ){
00056         muonSelect_= new StringCutObjectSelector<reco::Muon>(muonExtras.getParameter<std::string>("select"));
00057       }
00058       // isolation is optional; in case it's not found no
00059       // isolation will be applied
00060       if( muonExtras.existsAs<std::string>("isolation") ){
00061         muonIso_= new StringCutObjectSelector<reco::Muon>(muonExtras.getParameter<std::string>("isolation"));
00062       }
00063     }
00064     
00065     // jetExtras are optional; they may be omitted or 
00066     // empty
00067     if( cfg.existsAs<edm::ParameterSet>("jetExtras") ){
00068       edm::ParameterSet jetExtras=cfg.getParameter<edm::ParameterSet>("jetExtras");
00069       // jetCorrector is optional; in case it's not found 
00070       // the InputTag will remain empty
00071       if( jetExtras.existsAs<std::string>("jetCorrector") ){
00072         jetCorrector_= jetExtras.getParameter<std::string>("jetCorrector");
00073       }
00074       // read jetID information if it exists
00075       if(jetExtras.existsAs<edm::ParameterSet>("jetID")){
00076         edm::ParameterSet jetID=jetExtras.getParameter<edm::ParameterSet>("jetID");
00077         jetIDLabel_ =jetID.getParameter<edm::InputTag>("label");
00078         jetIDSelect_= new StringCutObjectSelector<reco::JetID>(jetID.getParameter<std::string>("select"));
00079       }
00080       // select is optional; in case it's not found no
00081       // selection will be applied (only implemented for 
00082       // CaloJets at the moment)
00083       if( jetExtras.existsAs<std::string>("select") ){
00084         jetSelect_= jetExtras.getParameter<std::string>("select");
00085       }
00086       // jetBDiscriminators are optional; in case they are
00087       // not found the InputTag will remain empty; they 
00088       // consist of pairs of edm::JetFlavorAssociation's & 
00089       // corresponding working points
00090       includeBTag_=jetExtras.existsAs<edm::ParameterSet>("jetBTaggers");
00091       if( includeBTag_ ){
00092         edm::ParameterSet btagEff=jetExtras.getParameter<edm::ParameterSet>("jetBTaggers").getParameter<edm::ParameterSet>("trackCountingEff");
00093         btagEff_= btagEff.getParameter<edm::InputTag>("label"); btagEffWP_= btagEff.getParameter<double>("workingPoint");
00094         edm::ParameterSet btagPur=jetExtras.getParameter<edm::ParameterSet>("jetBTaggers").getParameter<edm::ParameterSet>("trackCountingPur");
00095         btagPur_= btagPur.getParameter<edm::InputTag>("label"); btagPurWP_= btagPur.getParameter<double>("workingPoint");
00096         edm::ParameterSet btagVtx=jetExtras.getParameter<edm::ParameterSet>("jetBTaggers").getParameter<edm::ParameterSet>("secondaryVertex" );
00097         btagVtx_= btagVtx.getParameter<edm::InputTag>("label"); btagVtxWP_= btagVtx.getParameter<double>("workingPoint");
00098       }
00099     }
00100 
00101     // triggerExtras are optional; they may be omitted or empty
00102     if( cfg.existsAs<edm::ParameterSet>("triggerExtras") ){
00103       edm::ParameterSet triggerExtras=cfg.getParameter<edm::ParameterSet>("triggerExtras");
00104       triggerTable_=triggerExtras.getParameter<edm::InputTag>("src");
00105       triggerPaths_=triggerExtras.getParameter<std::vector<std::string> >("paths");
00106     }
00107 
00108     // massExtras is optional; in case it's not found no mass
00109     // window cuts are applied for the same flavor monitor
00110     // histograms
00111     if( cfg.existsAs<edm::ParameterSet>("massExtras") ){
00112       edm::ParameterSet massExtras=cfg.getParameter<edm::ParameterSet>("massExtras");
00113       lowerEdge_= massExtras.getParameter<double>("lowerEdge");
00114       upperEdge_= massExtras.getParameter<double>("upperEdge");
00115     }
00116 
00117     // setup the verbosity level for booking histograms;
00118     // per default the verbosity level will be set to 
00119     // STANDARD. This will also be the chosen level in
00120     // the case when the monitoring PSet is not found
00121     verbosity_=STANDARD;
00122     if( cfg.existsAs<edm::ParameterSet>("monitoring") ){
00123       edm::ParameterSet monitoring=cfg.getParameter<edm::ParameterSet>("monitoring");
00124       if(monitoring.getParameter<std::string>("verbosity") == "DEBUG"   )
00125         verbosity_= DEBUG;
00126       if(monitoring.getParameter<std::string>("verbosity") == "VERBOSE" )
00127         verbosity_= VERBOSE;
00128       if(monitoring.getParameter<std::string>("verbosity") == "STANDARD")
00129         verbosity_= STANDARD;
00130     }
00131     // and don't forget to do the histogram booking
00132     book(cfg.getParameter<std::string>("directory"));
00133   }
00134 
00135   void 
00136   MonitorEnsemble::book(std::string directory)
00137   {
00138     //set up the current directory path
00139     std::string current(directory); current+=label_;
00140     store_=edm::Service<DQMStore>().operator->();
00141     store_->setCurrentFolder(current);
00142 
00143     // determine number of bins for trigger monitoring
00144     unsigned int nPaths=triggerPaths_.size();
00145 
00146     // --- [STANDARD] --- //
00147     // pt of the leading muon
00148     hists_["muonPt_"     ] = store_->book1D("MuonPt"     , "pt(#mu)"          ,     50,     0.,    250.);   
00149     // muon multiplicity before std isolation
00150     hists_["muonMult_"   ] = store_->book1D("MuonMult"   , "N_{All}(#mu)"     ,     10,     0.,     10.);   
00151     // muon multiplicity after  std isolation
00152     hists_["muonMultIso_"] = store_->book1D("MuonMultIso", "N_{Iso}(#mu)"     ,     10,     0.,     10.);   
00153     // pt of the leading electron
00154     hists_["elecPt_"     ] = store_->book1D("ElecPt"     , "pt(e)"            ,     50,     0.,    250.);   
00155     // electron multiplicity before std isolation
00156     hists_["elecMult_"   ] = store_->book1D("ElecMult"   , "N_{All}(e)"       ,     10,     0.,     10.);   
00157     // electron multiplicity after  std isolation
00158     hists_["elecMultIso_"] = store_->book1D("ElecMultIso", "N_{Iso}(e)"       ,     10,     0.,     10.);   
00159     // multiplicity of jets with pt>20 (corrected to L2+L3)
00160     hists_["jetMult_"    ] = store_->book1D("JetMult"    , "N_{30}(jet)"      ,     10,     0.,     10.);   
00161     // trigger efficiency estimates for single lepton triggers
00162     hists_["triggerEff_" ] = store_->book1D("TriggerEff" , "Eff(trigger)"     , nPaths,     0.,  nPaths);
00163     // monitored trigger occupancy for single lepton triggers
00164     hists_["triggerMon_" ] = store_->book1D("TriggerMon" , "Mon(trigger)"     , nPaths,     0.,  nPaths);
00165     // MET (calo)
00166     hists_["metCalo_"    ] = store_->book1D("METCalo"    , "MET_{Calo}"       ,     50,     0.,    200.);   
00167     // W mass estimate
00168     hists_["massW_"      ] = store_->book1D("MassW"      , "M(W)"             ,     60,     0.,    300.);   
00169     // Top mass estimate
00170     hists_["massTop_"    ] = store_->book1D("MassTop"    , "M(Top)"           ,     50,     0.,    500.);   
00171 
00172     // set bin labels for trigger monitoring
00173     triggerBinLabels(std::string("trigger"), triggerPaths_);
00174 
00175     if( verbosity_==STANDARD) return;
00176 
00177     // --- [VERBOSE] --- //
00178     // eta of the leading muon
00179     hists_["muonEta_"    ] = store_->book1D("MuonEta"    , "#eta(#mu)"        ,     30,    -3.,      3.);   
00180     // std isolation variable of the leading muon
00181     hists_["muonRelIso_" ] = store_->book1D("MuonRelIso" , "Iso_{Rel}(#mu)"   ,     50,     0.,      1.);   
00182     // eta of the leading electron
00183     hists_["elecEta_"    ] = store_->book1D("ElecEta"    , "#eta(e)"          ,     30,    -3.,      3.);   
00184     // std isolation variable of the leading electron
00185     hists_["elecRelIso_" ] = store_->book1D("ElecRelIso" , "Iso_{Rel}(e)"     ,     50,     0.,      1.);   
00186     // multiplicity of btagged jets (for track counting high efficiency) with pt(L2L3)>30
00187     hists_["jetMultBEff_"] = store_->book1D("JetMultBEff", "N_{30}(b/eff)"    ,     10,     0.,     10.);   
00188     // btag discriminator for track counting high efficiency for jets with pt(L2L3)>30
00189     hists_["jetBDiscEff_"] = store_->book1D("JetBDiscEff", "Disc_{b/eff}(jet)",     100,     0.,     10.);   
00190     // pt of the 1. leading jet (corrected to L2+L3)
00191     hists_["jet1Pt_"     ] = store_->book1D("Jet1Pt"     , "pt_{L2L3}(jet1)"  ,     60,     0.,    300.);   
00192     // pt of the 2. leading jet (corrected to L2+L3)
00193     hists_["jet2Pt_"     ] = store_->book1D("Jet2Pt"     , "pt_{L2L3}(jet2)"  ,     60,     0.,    300.);   
00194     // pt of the 3. leading jet (corrected to L2+L3)
00195     hists_["jet3Pt_"     ] = store_->book1D("Jet3Pt"     , "pt_{L2L3}(jet3)"  ,     60,     0.,    300.);   
00196     // pt of the 4. leading jet (corrected to L2+L3)
00197     hists_["jet4Pt_"     ] = store_->book1D("Jet4Pt"     , "pt_{L2L3}(jet4)"  ,     60,     0.,    300.);   
00198     // MET (tc)
00199     hists_["metTC_"      ] = store_->book1D("METTC"      , "MET_{TC}"         ,     50,     0.,    200.);   
00200     // MET (pflow)
00201     hists_["metPflow_"   ] = store_->book1D("METPflow"   , "MET_{Pflow}"      ,     50,     0.,    200.);   
00202     // dz for muons (to suppress cosmis)
00203     hists_["muonDelZ_"    ] = store_->book1D("MuonDelZ"  , "d_{z}(#mu)"       ,     50,   -25.,     25.);
00204     // dxy for muons (to suppress cosmics)
00205     hists_["muonDelXY_"   ] = store_->book2D("MuonDelXY" , "d_{xy}(#mu)"      ,     50,   -0.1,     0.1,   50,   -0.1,   0.1);
00206 
00207     // set axes titles for dxy for muons
00208     hists_["muonDelXY_"   ]->setAxisTitle( "x [cm]", 1); hists_["muonDelXY_"   ]->setAxisTitle( "y [cm]", 2);
00209 
00210     if( verbosity_==VERBOSE) return;
00211 
00212     // --- [DEBUG] --- //
00213     // relative muon isolation in tracker for the leading muon
00214     hists_["muonTrkIso_" ] = store_->book1D("MuonTrkIso" , "Iso_{Trk}(#mu)"   ,     50,     0.,      1.);   
00215     // relative muon isolation in ecal+hcal for the leading muon
00216     hists_["muonCalIso_" ] = store_->book1D("MuonCalIso" , "Iso_{Ecal}(#mu)"  ,     50,     0.,      1.);   
00217     // relative electron isolation in tracker for the leading electron
00218     hists_["elecTrkIso_" ] = store_->book1D("ElecTrkIso" , "Iso_{Trk}(e)"     ,     50,     0.,      1.);   
00219     // relative electron isolation in ecal+hcal for the leading electron
00220     hists_["elecCalIso_" ] = store_->book1D("ElecCalIso" , "Iso_{Ecal}(e)"    ,     50,     0.,      1.);   
00221     // multiplicity of btagged jets (for track counting high purity) with pt(L2L3)>30
00222     hists_["jetMultBPur_"] = store_->book1D("JetMultBPur", "N_{30}(b/pur)"    ,     10,     0.,     10.);   
00223     // btag discriminator for track counting high purity
00224     hists_["jetBDiscPur_"] = store_->book1D("JetBDiscPur", "Disc_{b/pur}(Jet)",     100,     0.,     10.);   
00225     // multiplicity of btagged jets (for simple secondary vertex) with pt(L2L3)>30
00226     hists_["jetMultBVtx_"] = store_->book1D("JetMultBVtx", "N_{30}(b/vtx)"    ,     10,     0.,     10.);   
00227     // btag discriminator for simple secondary vertex
00228     hists_["jetBDiscVtx_"] = store_->book1D("JetBDiscVtx", "Disc_{b/vtx}(Jet)",     35,    -1.,      6.);   
00229     // pt of the 1. leading jet (uncorrected)
00230     hists_["jet1PtRaw_"  ] = store_->book1D("Jet1PtRaw"  , "pt_{Raw}(jet1)"   ,     60,     0.,    300.);   
00231     // pt of the 2. leading jet (uncorrected)
00232     hists_["jet2PtRaw_"  ] = store_->book1D("Jet2PtRaw"  , "pt_{Raw}(jet2)"   ,     60,     0.,    300.);   
00233     // pt of the 3. leading jet (uncorrected)
00234     hists_["jet3PtRaw_"  ] = store_->book1D("Jet3PtRaw"  , "pt_{Raw}(jet3)"   ,     60,     0.,    300.);   
00235     // pt of the 4. leading jet (uncorrected)
00236     hists_["jet4PtRaw_"  ] = store_->book1D("Jet4PtRaw"  , "pt_{Raw}(jet4)"   ,     60,     0.,    300.);   
00237     // selected events
00238     hists_["eventLogger_"] = store_->book2D("EventLogger", "Logged Events"    ,      9,     0.,      9.,   10,   0.,   10.);
00239 
00240     // set axes titles for selected events
00241     hists_["eventLogger_"]->getTH1()->SetOption("TEXT");
00242     hists_["eventLogger_"]->setBinLabel( 1 , "Run"             , 1);
00243     hists_["eventLogger_"]->setBinLabel( 2 , "Block"           , 1);
00244     hists_["eventLogger_"]->setBinLabel( 3 , "Event"           , 1);
00245     hists_["eventLogger_"]->setBinLabel( 4 , "pt_{L2L3}(jet1)" , 1);
00246     hists_["eventLogger_"]->setBinLabel( 5 , "pt_{L2L3}(jet2)" , 1);
00247     hists_["eventLogger_"]->setBinLabel( 6 , "pt_{L2L3}(jet3)" , 1);
00248     hists_["eventLogger_"]->setBinLabel( 7 , "pt_{L2L3}(jet4)" , 1);
00249     hists_["eventLogger_"]->setBinLabel( 8 , "M_{W}"           , 1);
00250     hists_["eventLogger_"]->setBinLabel( 9 , "M_{Top}"         , 1);
00251     hists_["eventLogger_"]->setAxisTitle("logged evts"         , 2);
00252     return;
00253   }
00254 
00255   void 
00256   MonitorEnsemble::fill(const edm::Event& event, const edm::EventSetup& setup)
00257   {
00258     // fetch trigger event if configured such 
00259     edm::Handle<edm::TriggerResults> triggerTable;
00260     if(!triggerTable_.label().empty()) {
00261       if( !event.getByLabel(triggerTable_, triggerTable) ) return;
00262     }
00263 
00264     /* 
00265     ------------------------------------------------------------
00266 
00267     Electron Monitoring
00268 
00269     ------------------------------------------------------------
00270     */
00271 
00272     // fill monitoring plots for electrons
00273     edm::Handle<edm::View<reco::GsfElectron> > elecs;
00274     if( !event.getByLabel(elecs_, elecs) ) return;
00275 
00276     // check availability of electron id
00277     edm::Handle<edm::ValueMap<float> > electronId; 
00278     if(!electronId_.label().empty()) {
00279       if( !event.getByLabel(electronId_, electronId) ) return;
00280     }
00281 
00282     // loop electron collection
00283     unsigned int eMult=0, eMultIso=0;
00284     std::vector<const reco::GsfElectron*> isoElecs;
00285     for(edm::View<reco::GsfElectron>::const_iterator elec=elecs->begin(); elec!=elecs->end(); ++elec){
00286       unsigned int idx = elec-elecs->begin();
00287       // restrict to electrons with good electronId
00288       if( electronId_.label().empty() ? true : ((int)(*electronId)[elecs->refAt(idx)] & eidPattern_) ){
00289         if(!elecSelect_ || (*elecSelect_)(*elec)){
00290           double isolationTrk = elec->pt()/(elec->pt()+elec->dr03TkSumPt());
00291           double isolationCal = elec->pt()/(elec->pt()+elec->dr03EcalRecHitSumEt()+elec->dr03HcalTowerSumEt());
00292           double isolationRel = (elec->dr03TkSumPt()+elec->dr03EcalRecHitSumEt()+elec->dr03HcalTowerSumEt())/elec->pt();
00293           if( eMult==0 ){
00294             // restrict to the leading electron
00295             fill("elecPt_" , elec->pt() );
00296             fill("elecEta_", elec->eta());
00297             fill("elecRelIso_" , isolationRel );
00298             fill("elecTrkIso_" , isolationTrk );
00299             fill("elecCalIso_" , isolationCal );
00300           }
00301           // in addition to the multiplicity counter buffer the iso 
00302           // electron candidates for later overlap check with jets
00303           ++eMult; if(!elecIso_ || (*elecIso_)(*elec)){ isoElecs.push_back(&(*elec)); ++eMultIso;}
00304         }
00305       }
00306     }
00307     fill("elecMult_",    eMult   );
00308     fill("elecMultIso_", eMultIso);
00309     
00310     /* 
00311     ------------------------------------------------------------
00312 
00313     Muon Monitoring
00314 
00315     ------------------------------------------------------------
00316     */
00317 
00318     // fill monitoring plots for muons
00319     unsigned int mMult=0, mMultIso=0;
00320 
00321     edm::Handle<edm::View<reco::Muon> > muons;
00322     if( !event.getByLabel(muons_, muons) ) return;
00323 
00324     for(edm::View<reco::Muon>::const_iterator muon=muons->begin(); muon!=muons->end(); ++muon){
00325       // restrict to globalMuons
00326       if( muon->isGlobalMuon() ){ 
00327         fill("muonDelZ_" , muon->globalTrack()->vz());
00328         fill("muonDelXY_", muon->globalTrack()->vx(), muon->globalTrack()->vy());
00329         // apply preselection
00330         if(!muonSelect_ || (*muonSelect_)(*muon)){
00331           double isolationTrk = muon->pt()/(muon->pt()+muon->isolationR03().sumPt);
00332           double isolationCal = muon->pt()/(muon->pt()+muon->isolationR03().emEt+muon->isolationR03().hadEt);
00333           double isolationRel = (muon->isolationR03().sumPt+muon->isolationR03().emEt+muon->isolationR03().hadEt)/muon->pt();
00334           if( mMult==0 ){
00335             // restrict to leading muon
00336             fill("muonPt_"     , muon->pt() );
00337             fill("muonEta_"    , muon->eta());
00338             fill("muonRelIso_" , isolationRel );
00339             fill("muonTrkIso_" , isolationTrk );
00340             fill("muonCalIso_" , isolationCal );
00341           }
00342            ++mMult; if(!muonIso_ || (*muonIso_)(*muon)) ++mMultIso;
00343         }
00344       }
00345     }
00346     fill("muonMult_",    mMult   );
00347     fill("muonMultIso_", mMultIso);
00348 
00349     /* 
00350     ------------------------------------------------------------
00351 
00352     Jet Monitoring
00353 
00354     ------------------------------------------------------------
00355     */
00356 
00357     // check availability of the btaggers
00358     edm::Handle<reco::JetTagCollection> btagEff, btagPur, btagVtx;
00359     if( includeBTag_ ){ 
00360       if( !event.getByLabel(btagEff_, btagEff) ) return;
00361       if( !event.getByLabel(btagPur_, btagPur) ) return;
00362       if( !event.getByLabel(btagVtx_, btagVtx) ) return;
00363     }
00364     // load jet corrector if configured such
00365     const JetCorrector* corrector=0;
00366     if(!jetCorrector_.empty()){
00367       // check whether a jet correcto is in the event setup or not
00368       if(setup.find( edm::eventsetup::EventSetupRecordKey::makeKey<JetCorrectionsRecord>() )){
00369         corrector = JetCorrector::getJetCorrector(jetCorrector_, setup);
00370       }
00371       else{ 
00372         edm::LogVerbatim( "TopSingleLeptonDQM" ) 
00373           << "\n"
00374           << "------------------------------------------------------------------------------------- \n"
00375           << " No JetCorrectionsRecord available from EventSetup:                                   \n" 
00376           << "  - Jets will not be corrected.                                                       \n"
00377           << "  - If you want to change this add the following lines to your cfg file:              \n"
00378           << "                                                                                      \n"
00379           << "  ## load jet corrections                                                             \n"
00380           << "  process.load(\"JetMETCorrections.Configuration.JetCorrectionServicesAllAlgos_cff\") \n"
00381           << "  process.prefer(\"ak5CaloL2L3\")                                                     \n"
00382           << "                                                                                      \n"
00383           << "------------------------------------------------------------------------------------- \n";
00384       }
00385     }
00386 
00387     // loop jet collection
00388     std::vector<reco::Jet> correctedJets;
00389     unsigned int mult=0, multBEff=0, multBPur=0, multBVtx=0;
00390 
00391     edm::Handle<edm::View<reco::Jet> > jets; 
00392     if( !event.getByLabel(jets_, jets) ) return;
00393 
00394     edm::Handle<reco::JetIDValueMap> jetID; 
00395     if(jetIDSelect_){ 
00396       if( !event.getByLabel(jetIDLabel_, jetID) ) return;
00397     }
00398 
00399     for(edm::View<reco::Jet>::const_iterator jet=jets->begin(); jet!=jets->end(); ++jet){
00400       // check jetID for calo jets
00401       unsigned int idx = jet-jets->begin();
00402       if( jetIDSelect_ && dynamic_cast<const reco::CaloJet*>(jets->refAt(idx).get())){
00403         if(!(*jetIDSelect_)((*jetID)[jets->refAt(idx)])) continue;
00404       }
00405       // chekc additional jet selection for calo, pf and bare reco jets
00406       if(dynamic_cast<const reco::CaloJet*>(&*jet)){
00407         reco::CaloJet sel = dynamic_cast<const reco::CaloJet&>(*jet); sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00408         StringCutObjectSelector<reco::CaloJet> jetSelect(jetSelect_); if(!jetSelect(sel)){ continue;}
00409       }
00410       else if(dynamic_cast<const reco::PFJet*>(&*jet)){
00411         reco::PFJet sel= dynamic_cast<const reco::PFJet&>(*jet); sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00412         StringCutObjectSelector<reco::PFJet> jetSelect(jetSelect_); if(!jetSelect(sel)) continue;
00413       } 
00414       else{
00415         reco::Jet sel = *jet; sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00416         StringCutObjectSelector<reco::Jet> jetSelect(jetSelect_); if(!jetSelect(sel)) continue;
00417       }
00418       // check for overlaps -- comment this to be synchronous with the selection
00419       //bool overlap=false;
00420       //for(std::vector<const reco::GsfElectron*>::const_iterator elec=isoElecs.begin(); elec!=isoElecs.end(); ++elec){
00421       //  if(reco::deltaR((*elec)->eta(), (*elec)->phi(), jet->eta(), jet->phi())<0.4){overlap=true; break;}
00422       //} if(overlap){continue;}
00423 
00424       // prepare jet to fill monitor histograms
00425       reco::Jet monitorJet = *jet; monitorJet.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00426       correctedJets.push_back(monitorJet);
00427       ++mult; // determine jet multiplicity
00428       if( includeBTag_ ){
00429         // fill b-discriminators
00430         edm::RefToBase<reco::Jet> jetRef = jets->refAt(idx);    
00431         fill("jetBDiscEff_", (*btagEff)[jetRef]); if( (*btagEff)[jetRef]>btagEffWP_ ) ++multBEff; 
00432         fill("jetBDiscPur_", (*btagPur)[jetRef]); if( (*btagPur)[jetRef]>btagPurWP_ ) ++multBPur; 
00433         fill("jetBDiscVtx_", (*btagVtx)[jetRef]); if( (*btagVtx)[jetRef]>btagVtxWP_ ) ++multBVtx; 
00434       }
00435       // fill pt (raw or L2L3) for the leading four jets  
00436       if(idx==0) {fill("jet1Pt_" , monitorJet.pt()); fill("jet1PtRaw_", jet->pt() );}
00437       if(idx==1) {fill("jet2Pt_" , monitorJet.pt()); fill("jet2PtRaw_", jet->pt() );}
00438       if(idx==2) {fill("jet3Pt_" , monitorJet.pt()); fill("jet3PtRaw_", jet->pt() );}
00439       if(idx==3) {fill("jet4Pt_" , monitorJet.pt()); fill("jet4PtRaw_", jet->pt() );}
00440     }
00441     fill("jetMult_"    , mult    );
00442     fill("jetMultBEff_", multBEff);
00443     fill("jetMultBPur_", multBPur);
00444     fill("jetMultBVtx_", multBVtx);
00445     
00446     /* 
00447     ------------------------------------------------------------
00448 
00449     MET Monitoring
00450 
00451     ------------------------------------------------------------
00452     */
00453 
00454     // fill monitoring histograms for met
00455     for(std::vector<edm::InputTag>::const_iterator met_=mets_.begin(); met_!=mets_.end(); ++met_){
00456       edm::Handle<edm::View<reco::MET> > met;
00457       if( !event.getByLabel(*met_, met) ) continue;
00458       if(met->begin()!=met->end()){
00459         unsigned int idx=met_-mets_.begin();
00460         if(idx==0) fill("metCalo_" , met->begin()->et());
00461         if(idx==1) fill("metTC_"   , met->begin()->et());
00462         if(idx==2) fill("metPflow_", met->begin()->et());
00463       }
00464     }
00465 
00466     /* 
00467     ------------------------------------------------------------
00468 
00469     Event Monitoring
00470 
00471     ------------------------------------------------------------
00472     */
00473 
00474     // fill W boson and top mass estimates
00475     Calculate eventKinematics(MAXJETS, WMASS);
00476     double wMass   = eventKinematics.massWBoson  (correctedJets);
00477     double topMass = eventKinematics.massTopQuark(correctedJets);
00478     if(wMass>=0 && topMass>=0) {fill("massW_" , wMass  ); fill("massTop_" , topMass);}
00479     // fill plots for trigger monitoring
00480     if((lowerEdge_==-1. && upperEdge_==-1.) || (lowerEdge_<wMass && wMass<upperEdge_) ){
00481       if(!triggerTable_.label().empty()) fill(event, *triggerTable, "trigger", triggerPaths_);
00482       if(logged_<=hists_.find("eventLogger_")->second->getNbinsY()){
00483         // log runnumber, lumi block, event number & some
00484         // more pysics infomation for interesting events
00485         fill("eventLogger_", 0.5, logged_+0.5, event.eventAuxiliary().run()); 
00486         fill("eventLogger_", 1.5, logged_+0.5, event.eventAuxiliary().luminosityBlock()); 
00487         fill("eventLogger_", 2.5, logged_+0.5, event.eventAuxiliary().event()); 
00488         if(correctedJets.size()>0) fill("eventLogger_", 3.5, logged_+0.5, correctedJets[0].pt()); 
00489         if(correctedJets.size()>1) fill("eventLogger_", 4.5, logged_+0.5, correctedJets[1].pt()); 
00490         if(correctedJets.size()>2) fill("eventLogger_", 5.5, logged_+0.5, correctedJets[2].pt()); 
00491         if(correctedJets.size()>3) fill("eventLogger_", 6.5, logged_+0.5, correctedJets[3].pt()); 
00492         fill("eventLogger_", 7.5, logged_+0.5, wMass  ); 
00493         fill("eventLogger_", 8.5, logged_+0.5, topMass); 
00494         ++logged_;
00495       }
00496     }
00497   }
00498   
00499 }
00500 
00501 
00502 TopSingleLeptonDQM::TopSingleLeptonDQM(const edm::ParameterSet& cfg): triggerTable_(""), beamspot_("") 
00503 {
00504   // configure preselection
00505   edm::ParameterSet presel=cfg.getParameter<edm::ParameterSet>("preselection");
00506   if( presel.existsAs<edm::ParameterSet>("trigger") ){
00507     edm::ParameterSet trigger=presel.getParameter<edm::ParameterSet>("trigger");
00508     triggerTable_=trigger.getParameter<edm::InputTag>("src");
00509     triggerPaths_=trigger.getParameter<std::vector<std::string> >("select");
00510   } 
00511   if( presel.existsAs<edm::ParameterSet>("vertex" ) ){
00512     edm::ParameterSet vertex=presel.getParameter<edm::ParameterSet>("vertex");
00513     vertex_= vertex.getParameter<edm::InputTag>("src");
00514     vertexSelect_= new StringCutObjectSelector<reco::Vertex>(vertex.getParameter<std::string>("select"));
00515   }
00516   if( presel.existsAs<edm::ParameterSet>("beamspot" ) ){
00517     edm::ParameterSet beamspot=presel.getParameter<edm::ParameterSet>("beamspot");
00518     beamspot_= beamspot.getParameter<edm::InputTag>("src");
00519     beamspotSelect_= new StringCutObjectSelector<reco::BeamSpot>(beamspot.getParameter<std::string>("select"));
00520   }
00521 
00522   // conifgure the selection
00523   std::vector<edm::ParameterSet> sel=cfg.getParameter<std::vector<edm::ParameterSet> >("selection");
00524   for(unsigned int i=0; i<sel.size(); ++i){
00525     selectionOrder_.push_back(sel.at(i).getParameter<std::string>("label"));
00526     selection_[selectionStep(selectionOrder_.back())] = std::make_pair(sel.at(i), new TopSingleLepton::MonitorEnsemble(selectionStep(selectionOrder_.back()).c_str(), cfg.getParameter<edm::ParameterSet>("setup")));
00527   }
00528 }
00529 
00530 void 
00531 TopSingleLeptonDQM::analyze(const edm::Event& event, const edm::EventSetup& setup)
00532 { 
00533   if(!triggerTable_.label().empty()){
00534     edm::Handle<edm::TriggerResults> triggerTable;
00535     if( !event.getByLabel(triggerTable_, triggerTable) ) return;
00536     if(!accept(event, *triggerTable, triggerPaths_)) return;
00537   }
00538   if(!beamspot_.label().empty()){
00539     edm::Handle<reco::BeamSpot> beamspot;
00540     if( !event.getByLabel(beamspot_, beamspot) ) return;
00541     if(!(*beamspotSelect_)(*beamspot)) return;
00542   }
00543   // apply selection steps
00544   unsigned int passed=0;
00545   for(std::vector<std::string>::const_iterator selIt=selectionOrder_.begin(); selIt!=selectionOrder_.end(); ++selIt){
00546     std::string key = selectionStep(*selIt), type = objectType(*selIt);
00547     if(selection_.find(key)!=selection_.end()){
00548       if(type=="empty"){
00549         selection_[key].second->fill(event, setup);
00550       }
00551       if(type=="muons"){
00552         SelectionStep<reco::Muon> step(selection_[key].first);
00553         if(step.select(event)){ ++passed;
00554           selection_[key].second->fill(event, setup);
00555         } else break;
00556       }
00557       if(type=="elecs"){
00558         SelectionStep<reco::GsfElectron> step(selection_[key].first);
00559         if(step.select(event)){ ++passed;
00560           selection_[key].second->fill(event, setup);
00561         } else break;
00562       }
00563       if(type=="jets" ){
00564         SelectionStep<reco::Jet> step(selection_[key].first);
00565         if(step.select(event, setup)){ ++passed;
00566           selection_[key].second->fill(event, setup);
00567         } else break;
00568       }
00569       if(type=="jets/pf" ){
00570         SelectionStep<reco::PFJet> step(selection_[key].first);
00571         if(step.select(event, setup)){ ++passed;
00572           selection_[key].second->fill(event, setup);
00573         } else break;
00574       }
00575       if(type=="jets/calo" ){
00576         SelectionStep<reco::CaloJet> step(selection_[key].first);
00577         if(step.select(event, setup)){ ++passed;
00578           selection_[key].second->fill(event, setup);
00579         } else break;
00580       }
00581       if(type=="met" ){
00582         SelectionStep<reco::MET> step(selection_[key].first);
00583         if(step.select(event)){ ++passed;
00584           selection_[key].second->fill(event, setup);
00585         } else break;
00586       }
00587     }
00588   }
00589 }
00590 
00591