CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_2_9/src/HLTriggerOffline/Top/src/TopHLTDiLeptonOfflineDQM.cc

Go to the documentation of this file.
00001 //#include <algorithm>
00002 #include "JetMETCorrections/Objects/interface/JetCorrectionsRecord.h"
00003 #include "DataFormats/JetReco/interface/CaloJet.h"
00004 #include "DataFormats/BTauReco/interface/JetTag.h"
00005 #include "HLTriggerOffline/Top/src/TopHLTDiLeptonOfflineDQM.h"
00006 #include "DataFormats/JetReco/interface/PFJet.h"
00007 #include "HLTriggerOffline/Top/interface/TopHLTDQMHelper.h"
00008 /*Originally from DQM/Physics by R. Wolf and J. Andrea*/
00009 
00010 namespace TopDiLeptonOffline {
00011 
00012   MonitorEnsemble::MonitorEnsemble(const char* label, const edm::ParameterSet& cfg) : 
00013    label_(label), eidPattern_(0), elecIso_(0), elecSelect_(0), muonIso_(0), muonSelect_(0), jetIDSelect_(0), 
00014    lowerEdge_(-1.), upperEdge_(-1.), elecMuLogged_(0), diMuonLogged_(0), diElecLogged_(0)
00015   {
00016     // sources have to be given; this PSet is not optional
00017     edm::ParameterSet sources=cfg.getParameter<edm::ParameterSet>("sources");
00018     muons_= sources.getParameter<edm::InputTag>("muons");
00019     elecs_= sources.getParameter<edm::InputTag>("elecs");
00020     jets_ = sources.getParameter<edm::InputTag>("jets" );
00021     mets_ = sources.getParameter<std::vector<edm::InputTag> >("mets" );
00022 
00023     // elecExtras are optional; they may be omitted or empty
00024     if( cfg.existsAs<edm::ParameterSet>("elecExtras") ){
00025       edm::ParameterSet elecExtras=cfg.getParameter<edm::ParameterSet>("elecExtras");
00026       // select is optional; in case it's not found no
00027       // selection will be applied
00028       if( elecExtras.existsAs<std::string>("select") ){
00029         elecSelect_= new StringCutObjectSelector<reco::GsfElectron>(elecExtras.getParameter<std::string>("select"));
00030       }
00031       // isolation is optional; in case it's not found no
00032       // isolation will be applied
00033       if( elecExtras.existsAs<std::string>("isolation") ){
00034         elecIso_= new StringCutObjectSelector<reco::GsfElectron>(elecExtras.getParameter<std::string>("isolation"));
00035       }
00036       // electronId is optional; in case it's not found the 
00037       // InputTag will remain empty
00038       if( elecExtras.existsAs<edm::ParameterSet>("electronId") ){
00039         edm::ParameterSet elecId=elecExtras.getParameter<edm::ParameterSet>("electronId");
00040         electronId_= elecId.getParameter<edm::InputTag>("src");
00041         eidPattern_= elecId.getParameter<int>("pattern");
00042       }
00043     }
00044     // muonExtras are optional; they may be omitted or empty
00045     if( cfg.existsAs<edm::ParameterSet>("muonExtras") ){
00046       edm::ParameterSet muonExtras=cfg.getParameter<edm::ParameterSet>("muonExtras");
00047       // select is optional; in case it's not found no
00048       // selection will be applied
00049       if( muonExtras.existsAs<std::string>("select") ){
00050         muonSelect_= new StringCutObjectSelector<reco::Muon>(muonExtras.getParameter<std::string>("select"));
00051       }
00052       // isolation is optional; in case it's not found no
00053       // isolation will be applied
00054       if( muonExtras.existsAs<std::string>("isolation") ){
00055         muonIso_= new StringCutObjectSelector<reco::Muon>(muonExtras.getParameter<std::string>("isolation"));
00056       }
00057     }
00058     // jetExtras are optional; they may be omitted or empty
00059     if( cfg.existsAs<edm::ParameterSet>("jetExtras") ){
00060       edm::ParameterSet jetExtras=cfg.getParameter<edm::ParameterSet>("jetExtras");
00061       // jetCorrector is optional; in case it's not found 
00062       // the InputTag will remain empty
00063       if( jetExtras.existsAs<std::string>("jetCorrector") ){
00064         jetCorrector_= jetExtras.getParameter<std::string>("jetCorrector");
00065       }
00066       // read jetID information if it exists
00067       if(jetExtras.existsAs<edm::ParameterSet>("jetID")){
00068         edm::ParameterSet jetID=jetExtras.getParameter<edm::ParameterSet>("jetID");
00069         jetIDLabel_ =jetID.getParameter<edm::InputTag>("label");
00070         jetIDSelect_= new StringCutObjectSelector<reco::JetID>(jetID.getParameter<std::string>("select"));
00071       }
00072       // select is optional; in case it's not found no
00073       // selection will be applied (only implemented for 
00074       // CaloJets at the moment)
00075       if( jetExtras.existsAs<std::string>("select") ){
00076         jetSelect_= jetExtras.getParameter<std::string>("select");
00077       }
00078     }
00079     // triggerExtras are optional; they may be omitted or empty
00080     if( cfg.existsAs<edm::ParameterSet>("triggerExtras") ){
00081       edm::ParameterSet triggerExtras=cfg.getParameter<edm::ParameterSet>("triggerExtras");
00082       triggerTable_=triggerExtras.getParameter<edm::InputTag>("src");
00083       elecMuPaths_ =triggerExtras.getParameter<std::vector<std::string> >("pathsELECMU");
00084       diMuonPaths_ =triggerExtras.getParameter<std::vector<std::string> >("pathsDIMUON");
00085     }
00086     // massExtras is optional; in case it's not found no mass
00087     // window cuts are applied for the same flavor monitor
00088     // histograms
00089     if( cfg.existsAs<edm::ParameterSet>("massExtras") ){
00090       edm::ParameterSet massExtras=cfg.getParameter<edm::ParameterSet>("massExtras");
00091       lowerEdge_= massExtras.getParameter<double>("lowerEdge");
00092       upperEdge_= massExtras.getParameter<double>("upperEdge");
00093     }
00094 
00095     // setup the verbosity level for booking histograms;
00096     // per default the verbosity level will be set to 
00097     // STANDARD. This will also be the chosen level in
00098     // the case when the monitoring PSet is not found
00099     verbosity_=STANDARD;
00100     if( cfg.existsAs<edm::ParameterSet>("monitoring") ){
00101       edm::ParameterSet monitoring=cfg.getParameter<edm::ParameterSet>("monitoring");
00102       if(monitoring.getParameter<std::string>("verbosity") == "DEBUG"   )
00103         verbosity_= DEBUG;
00104       if(monitoring.getParameter<std::string>("verbosity") == "VERBOSE" )
00105         verbosity_= VERBOSE;
00106       if(monitoring.getParameter<std::string>("verbosity") == "STANDARD")
00107         verbosity_= STANDARD;
00108     }
00109     // and don't forget to do the histogram booking
00110     book(cfg.getParameter<std::string>("directory"));
00111   }
00112 
00113   void 
00114   MonitorEnsemble::book(std::string directory)
00115   {
00116     //set up the current directory path
00117     std::string current(directory); current+=label_;
00118     store_=edm::Service<DQMStore>().operator->();
00119     store_->setCurrentFolder(current);
00120 
00121     // determine number of bins for trigger monitoring
00122     unsigned int nElecMu=elecMuPaths_.size();
00123     unsigned int nDiMuon=diMuonPaths_.size();
00124 
00125     // --- [STANDARD] --- //
00126     //Run Number
00127     hists_["RunNumb_"    ] = store_->book1D("RunNumber"  , "Run Nr."          ,   1.e4,  1.5e5,    3.e5);
00128     // invariant mass of opposite charge lepton pair (only filled for same flavor)
00129     hists_["invMass_"     ] = store_->book1D("InvMass"     , "M(lep1, lep2)"           ,       80,   0.,     320.);
00130     // invariant mass of opposite charge lepton pair (only filled for same flavor)
00131     hists_["invMassLog_"  ] = store_->book1D("InvMassLog"  , "log_{10}(M(lep1, lep2))" ,       80,   .1,      2.5);
00132     // invariant mass of same charge lepton pair (log10 for low mass region, only filled for same flavor)
00133     hists_["invMassWC_"   ] = store_->book1D("InvMassWC"   , "M_{WC}(L1, L2)"          ,       80,   0.,     320.);
00134     // invariant mass of same charge lepton pair (log10 for low mass region, only filled for same flavor)
00135     hists_["invMassWCLog_"] = store_->book1D("InvMassLogWC", "log_{10}(M_{WC})"        ,       80,   .1,      2.5);
00136     // decay channel [1]: muon/muon, [2]:elec/elec, [3]:elec/muon 
00137     hists_["decayChannel_"] = store_->book1D("DecayChannel", "Decay Channel"           ,        3,    0,        3);
00138     // trigger efficiency estimates for the electron muon channel
00139     hists_["elecMuEff_"   ] = store_->book1D("ElecMuEff"   , "Eff(e/#mu paths)"        ,  nElecMu,   0.,  nElecMu);
00140     // monitored trigger occupancy for the electron muon channel
00141     hists_["elecMuMon_"   ] = store_->book1D("ElecMuMon"   , "Mon(e/#mu paths)"        ,  nElecMu,   0.,  nElecMu);
00142     // trigger efficiency estimates for the di muon channel
00143     hists_["diMuonEff_"   ] = store_->book1D("DiMuonEff"   , "Eff(#mu/#mu paths)"      ,  nDiMuon,   0.,  nDiMuon);
00144     // monitored trigger occupancy for the di muon channel
00145     hists_["diMuonMon_"   ] = store_->book1D("DiMuonMon"   , "Mon(#mu/#mu paths)"      ,  nDiMuon,   0.,  nDiMuon);
00146     // pt of the leading lepton
00147     hists_["lep1Pt_"      ] = store_->book1D("Lep1Pt"      , "pt(lep1)"                ,       50,   0.,     200.);
00148     // pt of the 2. leading lepton
00149     hists_["lep2Pt_"      ] = store_->book1D("Lep2Pt"      , "pt(lep2)"                ,       50,   0.,     200.);
00150     // multiplicity of jets with pt>30 (corrected to L2+L3)
00151     hists_["jetMult_"     ] = store_->book1D("JetMult"     , "N_{30}(jet)"             ,       21, -0.5,      20.5); 
00152     // MET (calo)
00153     hists_["metCalo_"     ] = store_->book1D("METCalo"     , "MET_{Calo}"              ,       50,   0.,     200.);
00154 
00155     // set bin labels for trigger monitoring
00156     triggerBinLabels(std::string("elecMu"), elecMuPaths_);
00157     triggerBinLabels(std::string("diMuon"), diMuonPaths_);
00158     // set bin labels for decayChannel_
00159     hists_["decayChannel_"]->setBinLabel( 1, "#mu e"  , 1);
00160     hists_["decayChannel_"]->setBinLabel( 2, "#mu #mu", 1);
00161     hists_["decayChannel_"]->setBinLabel( 3, "e e"    , 1);
00162 
00163     if( verbosity_==STANDARD) return;
00164 
00165     // --- [VERBOSE] --- //
00166     // mean eta of the candidate leptons
00167     hists_["sumEtaL1L2_"  ] = store_->book1D("SumEtaL1L2"  , "<#eta>(lep1, lep2)"      ,       100,  -5.,       5.); 
00168     // deltaEta between the 2 candidate leptons
00169     hists_["dEtaL1L2_"    ] = store_->book1D("DEtaL1L2"    , "#Delta#eta(lep1,lep2)"   ,       80,  -4.,       4.);
00170     // deltaPhi between the 2 candidate leptons
00171     hists_["dPhiL1L2_"    ] = store_->book1D("DPhiL1L2"    , "#Delta#phi(lep1,lep2)"   ,       64, -3.2,      3.2);
00172     // pt of the candidate electron (depending on the decay channel)
00173     hists_["elecPt_"      ] = store_->book1D("ElecPt"      , "pt(e)"                   ,       50,   0.,     200.);
00174     // relative isolation of the candidate electron (depending on the decay channel)
00175     hists_["elecRelIso_"  ] = store_->book1D("ElecRelIso"  , "Iso_{Rel}(e)"            ,       50,   0.,       1.);
00176     // pt of the canddiate muon (depending on the decay channel)
00177     hists_["muonPt_"      ] = store_->book1D("MuonPt"      , "pt(#mu)"                 ,       50,   0.,     200.);
00178     // relative isolation of the candidate muon (depending on the decay channel)
00179     hists_["muonRelIso_"  ] = store_->book1D("MuonRelIso"  , "Iso_{Rel}(#mu)"          ,       50,   0.,       1.);
00180     // pt of the 1. leading jet (corrected to L2+L3)
00181     hists_["jet1Pt_"      ] = store_->book1D("Jet1Pt"      , "pt_{L2L3}(jet1)"         ,       60,   0.,     300.);   
00182     // pt of the 2. leading jet (corrected to L2+L3)
00183     hists_["jet2Pt_"      ] = store_->book1D("Jet2Pt"      , "pt_{L2L3}(jet2)"         ,       60,   0.,     300.);
00184     // MET (PF)
00185     hists_["metPflow_"    ] = store_->book1D("METPflow"    , "MET_{Pflow}"             ,       50,   0.,     200.);
00186     // MET (TC)
00187     hists_["metTC_"       ] = store_->book1D("METTC"       , "MET_{TC}"                ,       50,   0.,     200.);
00188     // dz for muons (to suppress cosmis)
00189     hists_["muonDelZ_"    ] = store_->book1D("MuonDelZ"    , "d_{z}(#mu)"              ,       50, -25.,      25.);
00190     // dxy for muons (to suppress cosmics)
00191     hists_["muonDelXY_"   ] = store_->book2D("MuonDelXY"   , "d_{xy}(#mu)"             , 50,  -1.,  1., 50, -1., 1.);
00192     // lepton multiplicity after std isolation
00193     hists_["lepMultIso_"  ] = store_->book2D("LepMultIso"  , "N_{Iso}(e) vs N_{Iso}(#mu)" ,  5,    0.,   5.,  5,   0.,    5.);
00194 
00195     // set axes titles for dxy for muons
00196     hists_["muonDelXY_"   ]->setAxisTitle( "x [cm]", 1); hists_["muonDelXY_"   ]->setAxisTitle( "y [cm]", 2);
00197     // set axes titles for lepton multiplicity after std isolation
00198     hists_["lepMultIso_"  ]->setAxisTitle( "N_{Iso}(#mu)", 1); hists_["lepMultIso_"   ]->setAxisTitle( "N_{Iso}(elec)", 2);
00199 
00200     if( verbosity_==VERBOSE) return;
00201 
00202     // --- [DEBUG] --- //
00203     // electron multiplicity after std isolation
00204     hists_["elecMultIso_" ] = store_->book1D("ElecMultIso" , "N_{Iso}(e)"              ,       11, -0.5,     10.5);
00205     // muon multiplicity after std isolation
00206     hists_["muonMultIso_" ] = store_->book1D("MuonMultIso" , "N_{Iso}(#mu)"            ,       11, -0.5,     10.5);
00207     // calo isolation of the candidate muon (depending on the decay channel)
00208     hists_["muonCalIso_"  ] = store_->book1D("MuonCalIso"  , "Iso_{Cal}(#mu)"          ,       50,   0.,       1.);
00209     // track isolation of the candidate muon (depending on the decay channel)
00210     hists_["muonTrkIso_"  ] = store_->book1D("MuonTrkIso"  , "Iso_{Trk}(#mu)"          ,       50,   0.,       1.);
00211     // calo isolation of the candidate electron (depending on the decay channel)
00212     hists_["elecCalIso_"  ] = store_->book1D("ElecCalIso"  , "Iso_{Cal}(e)"            ,       50,   0.,       1.);
00213     // track isolation of the candidate electron (depending on the decay channel)
00214     hists_["elecTrkIso_"  ] = store_->book1D("ElecTrkIso"  , "Iso_{Trk}(e)"            ,       50,   0.,       1.);
00215     // eta of the leading jet
00216     hists_["jet1Eta_"     ] = store_->book1D("Jet1Eta"     , "#eta(jet1)"              ,       30,  -5.,       5.); 
00217     // eta of the 2. leading jet
00218     hists_["jet2Eta_"     ] = store_->book1D("Jet2Eta"     , "#eta(jet2)"              ,       30,  -5.,       5.);
00219     // pt of the 1. leading jet (not corrected)
00220     hists_["jet1PtRaw_"   ] = store_->book1D("Jet1PtRaw"   , "pt_{Raw}(jet1)"          ,       60,   0.,     300.);   
00221     // pt of the 2. leading jet (not corrected)     
00222     hists_["jet2PtRaw_"   ] = store_->book1D("Jet2PtRaw"   , "pt_{Raw}(jet2)"          ,       60,   0.,     300.);
00223     // deltaEta between the 2 leading jets
00224     hists_["dEtaJet1Jet2_"] = store_->book1D("DEtaJet1Jet2", "#Delta#eta(jet1,jet2)"   ,       80,  -4.,       4.);
00225     // deltaEta between the lepton and the leading jet
00226     hists_["dEtaJet1Lep1_"] = store_->book1D("DEtaJet1Lep1", "#Delta#eta(jet1,lep1)"   ,       80,  -4.,       4.);
00227     // deltaEta between the lepton and MET
00228     hists_["dEtaLep1MET_" ] = store_->book1D("DEtaLep1MET" , "#Delta#eta(lep1,MET)"    ,       80,  -4.,       4.);
00229     // deltaEta between leading jet and MET
00230     hists_["dEtaJet1MET_" ] = store_->book1D("DEtaJet1MET" , "#Delta#eta(jet1,MET)"    ,       80,  -4.,       4.);
00231     // deltaPhi of 2 leading jets
00232     hists_["dPhiJet1Jet2_"] = store_->book1D("DPhiJet1Jet2", "#Delta#phi(jet1,jet2)"   ,       64, -3.2,      3.2);
00233     // deltaPhi of 1. lepton and 1. jet
00234     hists_["dPhiJet1Lep1_"] = store_->book1D("DPhiJet1Lep1", "#Delta#phi(jet1,lep1)"   ,       64, -3.2,      3.2);
00235     // deltaPhi of 1. lepton and MET
00236     hists_["dPhiLep1MET_" ] = store_->book1D("DPhiLep1MET" , "#Delta#phi(lep1,MET)"    ,       64, -3.2,      3.2);
00237     // deltaPhi of 1. jet and MET
00238     hists_["dPhiJet1MET_" ] = store_->book1D("DPhiJet1MET" , "#Delta#phi(jet1,MET)"    ,       64, -3.2,      3.2);
00239     // selected dimuon events
00240     hists_["diMuonLogger_"] = store_->book2D("DiMuonLogger", "Logged DiMuon Events"    ,        8,   0.,       8.,   10,   0.,   10.);
00241     // selected dielec events
00242     hists_["diElecLogger_"] = store_->book2D("DiElecLogger", "Logged DiElec Events"    ,        8,   0.,       8.,   10,   0.,   10.);
00243     // selected elemu events
00244     hists_["elecMuLogger_"] = store_->book2D("ElecMuLogger", "Logged ElecMu Events"    ,        8,   0.,       8.,   10,   0.,   10.);
00245 
00246     // set bin labels for trigger monitoring
00247     loggerBinLabels(std::string("diMuonLogger_")); 
00248     loggerBinLabels(std::string("diElecLogger_")); 
00249     loggerBinLabels(std::string("elecMuLogger_"));
00250     return;
00251   }
00252 
00253   void 
00254   MonitorEnsemble::fill(const edm::Event& event, const edm::EventSetup& setup)
00255   {
00256     // fetch trigger event if configured such 
00257     edm::Handle<edm::TriggerResults> triggerTable;
00258     if(!triggerTable_.label().empty()) {
00259       if( !event.getByLabel(triggerTable_, triggerTable) ) return;
00260     }
00261 
00262     /*
00263     ------------------------------------------------------------
00264     
00265     Run and Inst. Luminosity information (Inst. Lumi. filled now with a dummy value=5.0)
00266     
00267     ------------------------------------------------------------
00268     */
00269     
00270     if (!event.eventAuxiliary().run()) return;
00271     fill("RunNumb_", event.eventAuxiliary().run());   
00272     
00273     double dummy=5.; fill("InstLumi_", dummy);
00274      
00275     
00276     /* 
00277     ------------------------------------------------------------
00278 
00279     Muon Selection
00280 
00281     ------------------------------------------------------------
00282     */
00283 
00284     // buffer isolated muons
00285     std::vector<const reco::Muon*> isoMuons;
00286 
00287     edm::Handle<edm::View<reco::Muon> > muons;
00288     if( !event.getByLabel(muons_, muons) ) return;
00289 
00290     for(edm::View<reco::Muon>::const_iterator muon=muons->begin(); muon!=muons->end(); ++muon){
00291       // restrict to globalMuons
00292       if( muon->isGlobalMuon() ){ 
00293         fill("muonDelZ_" , muon->globalTrack()->vz());
00294         fill("muonDelXY_", muon->globalTrack()->vx(), muon->globalTrack()->vy());
00295         // apply preselection
00296         if(!muonSelect_ || (*muonSelect_)(*muon)){
00297           double isolationTrk = muon->pt()/(muon->pt()+muon->isolationR03().sumPt);
00298           double isolationCal = muon->pt()/(muon->pt()+muon->isolationR03().emEt+muon->isolationR03().hadEt);
00299           double isolationRel = (muon->isolationR03().sumPt+muon->isolationR03().emEt+muon->isolationR03().hadEt)/muon->pt();
00300           fill("muonTrkIso_" , isolationTrk); fill("muonCalIso_" , isolationCal); fill("muonRelIso_" , isolationRel);
00301           if(!muonIso_ || (*muonIso_)(*muon)) isoMuons.push_back(&(*muon));
00302         }
00303       }
00304     }
00305     fill("muonMultIso_", isoMuons.size());
00306 
00307     /* 
00308     ------------------------------------------------------------
00309 
00310     Electron Selection
00311 
00312     ------------------------------------------------------------
00313     */
00314 
00315     // buffer isolated electronss
00316     std::vector<const reco::GsfElectron*> isoElecs;
00317     edm::Handle<edm::ValueMap<float> > electronId; 
00318     if(!electronId_.label().empty()) {
00319       if( !event.getByLabel(electronId_, electronId) ) return;
00320     }
00321 
00322     edm::Handle<edm::View<reco::GsfElectron> > elecs;
00323     if( !event.getByLabel(elecs_, elecs) ) return;
00324 
00325     for(edm::View<reco::GsfElectron>::const_iterator elec=elecs->begin(); elec!=elecs->end(); ++elec){
00326       // restrict to electrons with good electronId
00327       int idx = elec-elecs->begin();
00328       if( electronId_.label().empty() ? true : ((int)(*electronId)[elecs->refAt(idx)] & eidPattern_) ){
00329         // apply preselection
00330         if(!elecSelect_ || (*elecSelect_)(*elec)){
00331           double isolationTrk = elec->pt()/(elec->pt()+elec->dr03TkSumPt());
00332           double isolationCal = elec->pt()/(elec->pt()+elec->dr03EcalRecHitSumEt()+elec->dr03HcalTowerSumEt());
00333           double isolationRel = (elec->dr03TkSumPt()+elec->dr03EcalRecHitSumEt()+elec->dr03HcalTowerSumEt())/elec->pt();
00334           fill("elecTrkIso_" , isolationTrk); fill("elecCalIso_" , isolationCal); fill("elecRelIso_" , isolationRel);
00335           if(!elecIso_ || (*elecIso_)(*elec)) isoElecs.push_back(&(*elec));
00336         }
00337       }
00338     }
00339     fill("elecMultIso_", isoElecs.size());
00340 
00341     /* 
00342     ------------------------------------------------------------
00343 
00344     Jet Selection
00345 
00346     ------------------------------------------------------------
00347     */
00348 
00349     const JetCorrector* corrector=0;
00350     if(!jetCorrector_.empty()){
00351       // check whether a jet correcto is in the event setup or not
00352       if(setup.find( edm::eventsetup::EventSetupRecordKey::makeKey<JetCorrectionsRecord>() )){
00353         corrector = JetCorrector::getJetCorrector(jetCorrector_, setup);
00354       }
00355       else{
00356         edm::LogVerbatim( "TopHLTDiLeptonOfflineDQM" ) 
00357           << "\n"
00358           << "------------------------------------------------------------------------------------- \n"
00359           << " No JetCorrectionsRecord available from EventSetup:                                   \n" 
00360           << "  - Jets will not be corrected.                                                       \n"
00361           << "  - If you want to change this add the following lines to your cfg file:              \n"
00362           << "                                                                                      \n"
00363           << "  ## load jet corrections                                                             \n"
00364           << "  process.load(\"JetMETCorrections.Configuration.JetCorrectionServicesAllAlgos_cff\") \n"
00365           << "  process.prefer(\"ak5CaloL2L3\")                                                     \n"
00366           << "                                                                                      \n"
00367           << "------------------------------------------------------------------------------------- \n";
00368       }
00369     }
00370 
00371     unsigned int mult=0;
00372     // buffer leadingJets
00373     std::vector<reco::Jet> leadingJets;
00374     edm::Handle<edm::View<reco::Jet> > jets; 
00375     if( !event.getByLabel(jets_, jets) ) return;
00376 
00377     edm::Handle<reco::JetIDValueMap> jetID;
00378     if(jetIDSelect_){ 
00379       if( !event.getByLabel(jetIDLabel_, jetID) ) return;
00380     }
00381 
00382     for(edm::View<reco::Jet>::const_iterator jet=jets->begin(); jet!=jets->end(); ++jet){
00383       unsigned int idx=jet-jets->begin();
00384       if( jetIDSelect_ && dynamic_cast<const reco::CaloJet*>(jets->refAt(idx).get())){
00385         if(!(*jetIDSelect_)((*jetID)[jets->refAt(idx)])) continue;
00386       }
00387       // chekc additional jet selection for calo, pf and bare reco jets
00388       if(dynamic_cast<const reco::CaloJet*>(&*jet)){
00389         reco::CaloJet sel = dynamic_cast<const reco::CaloJet&>(*jet); sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00390         StringCutObjectSelector<reco::CaloJet> jetSelect(jetSelect_); if(!jetSelect(sel)){ continue;}
00391       }
00392       else if(dynamic_cast<const reco::PFJet*>(&*jet)){
00393         reco::PFJet sel= dynamic_cast<const reco::PFJet&>(*jet); sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00394         StringCutObjectSelector<reco::PFJet> jetSelect(jetSelect_); if(!jetSelect(sel)) continue;
00395       } 
00396       else{
00397         reco::Jet sel = *jet; sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
00398         StringCutObjectSelector<reco::Jet> jetSelect(jetSelect_); if(!jetSelect(sel)) continue;
00399       }
00400       // check for overlaps
00401       bool overlap=false;
00402       for(std::vector<const reco::GsfElectron*>::const_iterator elec=isoElecs.begin(); elec!=isoElecs.end(); ++elec){
00403         if(reco::deltaR((*elec)->eta(), (*elec)->phi(), jet->eta(), jet->phi())<0.4){overlap=true; break;}
00404       } if(overlap){continue;}
00405       // prepare jet to fill monitor histograms
00406       reco::Jet monitorJet=*jet; monitorJet.scaleEnergy(corrector ?  corrector->correction(*jet) : 1.);
00407       ++mult; // determine jet multiplicity
00408       if(idx==0) {
00409         leadingJets.push_back(monitorJet);
00410         fill("jet1Pt_"    , monitorJet.pt());
00411         fill("jet1PtRaw_" , jet->pt() );
00412         fill("jet1Eta_"   , jet->eta());
00413       }
00414       if(idx==1) {
00415         leadingJets.push_back(monitorJet);
00416         fill("jet2Pt_"    , monitorJet.pt());
00417         fill("jet2PtRaw_" , jet->pt() );
00418         fill("jet2Eta_"   , jet->eta());
00419       }
00420     }
00421     if(leadingJets.size()>1){
00422       fill("dEtaJet1Jet2_" , leadingJets[0].eta()-leadingJets[1].eta());
00423       fill("dPhiJet1Jet2_" , reco::deltaPhi(leadingJets[0].phi(), leadingJets[1].phi()));
00424       if( !isoMuons.empty() ){
00425         if( isoElecs.empty() || isoMuons[0]->pt()>isoElecs[0]->pt() ){
00426           fill("dEtaJet1Lep1_" , isoMuons[0]->eta()-leadingJets[0].eta());
00427           fill("dPhiJet1Lep1_" , reco::deltaPhi(isoMuons[0]->phi() , leadingJets[0].phi()));
00428         } 
00429       }
00430       if( !isoElecs.empty() ){
00431         if( isoMuons.empty() || isoElecs[0]->pt()>isoMuons[0]->pt() ){
00432           fill("dEtaJet1Lep1_" , isoElecs[0]->eta()-leadingJets[0].eta());
00433           fill("dPhiJet1Lep1_" , reco::deltaPhi(isoElecs[0]->phi() , leadingJets[0].phi()));
00434         }
00435       }
00436     }
00437     fill("jetMult_", mult);
00438     
00439     /* 
00440     ------------------------------------------------------------
00441 
00442     MET Selection
00443 
00444     ------------------------------------------------------------
00445     */
00446 
00447     // buffer for event logging 
00448     reco::MET caloMET;
00449     for(std::vector<edm::InputTag>::const_iterator met_=mets_.begin(); met_!=mets_.end(); ++met_){
00450 
00451       edm::Handle<edm::View<reco::MET> > met;
00452       if( !event.getByLabel(*met_, met) ) continue;
00453 
00454       if(met->begin()!=met->end()){
00455         unsigned int idx=met_-mets_.begin();
00456         if(idx==0){
00457           caloMET=*met->begin(); 
00458           fill("metCalo_", met->begin()->et());
00459           if(!leadingJets.empty()){
00460             fill("dEtaJet1MET_" , leadingJets[0].eta()-met->begin()->eta());
00461             fill("dPhiJet1MET_" , reco::deltaPhi(leadingJets[0].phi(), met->begin()->phi()));
00462           }
00463           if( !isoMuons.empty() ){
00464             if( isoElecs.empty() || isoMuons[0]->pt()>isoElecs[0]->pt() ){
00465               fill("dEtaLep1MET_" , isoMuons[0]->eta()-met->begin()->eta());
00466               fill("dPhiLep1MET_" , reco::deltaPhi(isoMuons[0]->phi(), met->begin()->phi()));
00467             } 
00468           }
00469           if( !isoElecs.empty() ){
00470             if( isoMuons.empty() || isoElecs[0]->pt()>isoMuons[0]->pt() ){
00471               fill("dEtaLep1MET_" , isoElecs[0]->eta()-met->begin()->eta());
00472               fill("dPhiLep1MET_" , reco::deltaPhi(isoElecs[0]->phi(), met->begin()->phi()));
00473             }
00474           }
00475         }
00476         if(idx==1){ fill("metTC_"   , met->begin()->et());}
00477         if(idx==2){ fill("metPflow_", met->begin()->et());}
00478       }
00479     }
00480 
00481 
00482     /* 
00483     ------------------------------------------------------------
00484 
00485     Event Monitoring
00486 
00487     ------------------------------------------------------------
00488     */
00489 
00490     // check number of isolated leptons
00491     fill("lepMultIso_", isoMuons.size(), isoElecs.size());
00492     // ELECMU channel
00493     if( decayChannel(isoMuons, isoElecs) == ELECMU ){
00494       fill("decayChannel_", 0.5);
00495       double mass = (isoElecs[0]->p4()+isoMuons[0]->p4()).mass();
00496       if( (lowerEdge_==-1. && upperEdge_==-1.) || (lowerEdge_<mass && mass<upperEdge_) ){
00497         
00498         fill("dEtaL1L2_"  , isoElecs[0]->eta()-isoMuons[0]->eta()); 
00499         fill("sumEtaL1L2_", (isoElecs[0]->eta()+isoMuons[0]->eta())/2); 
00500         fill("dPhiL1L2_"  , reco::deltaPhi(isoElecs[0]->phi(), isoMuons[0]->eta())); 
00501         fill("elecPt_", isoElecs[0]->pt()); fill("muonPt_", isoMuons[0]->pt()); 
00502         fill("lep1Pt_", isoElecs[0]->pt()>isoMuons[0]->pt() ? isoElecs[0]->pt() : isoMuons[0]->pt());
00503         fill("lep2Pt_", isoElecs[0]->pt()>isoMuons[0]->pt() ? isoMuons[0]->pt() : isoElecs[0]->pt());
00504         // fill plots for trigger monitoring
00505         if(!triggerTable_.label().empty()) fill(event, *triggerTable, "elecMu", elecMuPaths_);
00506         if(elecMuLogged_<=hists_.find("elecMuLogger_")->second->getNbinsY()){
00507           // log runnumber, lumi block, event number & some
00508           // more pysics infomation for interesting events
00509           fill("elecMuLogger_", 0.5, elecMuLogged_+0.5, event.eventAuxiliary().run()); 
00510           fill("elecMuLogger_", 1.5, elecMuLogged_+0.5, event.eventAuxiliary().luminosityBlock()); 
00511           fill("elecMuLogger_", 2.5, elecMuLogged_+0.5, event.eventAuxiliary().event()); 
00512           fill("elecMuLogger_", 3.5, elecMuLogged_+0.5, isoMuons[0]->pt()); 
00513           fill("elecMuLogger_", 4.5, elecMuLogged_+0.5, isoElecs[0]->pt()); 
00514           if(leadingJets.size()>0) fill("elecMuLogger_", 5.5, elecMuLogged_+0.5, leadingJets[0].pt()); 
00515           if(leadingJets.size()>1) fill("elecMuLogger_", 6.5, elecMuLogged_+0.5, leadingJets[1].pt()); 
00516           fill("elecMuLogger_", 7.5, elecMuLogged_+0.5, caloMET.et()); 
00517           ++elecMuLogged_; 
00518         }
00519       }
00520     }
00521 
00522     // DIMUON channel
00523     if( decayChannel(isoMuons, isoElecs) == DIMUON ){
00524       fill("decayChannel_", 1.5);
00525       int charge = isoMuons[0]->charge()*isoMuons[1]->charge();
00526       double mass = (isoMuons[0]->p4()+isoMuons[1]->p4()).mass();
00527       
00528       fill(charge<0 ? "invMass_"    : "invMassWC_"    , mass       );
00529       fill(charge<0 ? "invMassLog_" : "invMassWCLog_" , log10(mass));
00530       if((lowerEdge_==-1. && upperEdge_==-1.) || (lowerEdge_<mass && mass<upperEdge_) ){
00531         fill("dEtaL1L2_"  , isoMuons[0]->eta()-isoMuons[1]->eta() );
00532         fill("sumEtaL1L2_", (isoMuons[0]->eta()+isoMuons[1]->eta())/2);
00533         fill("dPhiL1L2_", reco::deltaPhi(isoMuons[0]->phi(),isoMuons[1]->phi()) );
00534         fill("muonPt_", isoMuons[0]->pt()); fill("muonPt_", isoMuons[1]->pt()); 
00535         fill("lep1Pt_", isoMuons[0]->pt()); fill("lep2Pt_", isoMuons[1]->pt()); 
00536         // fill plots for trigger monitoring
00537         if(!triggerTable_.label().empty()) fill(event, *triggerTable, "diMuon", diMuonPaths_);
00538         if(diMuonLogged_<=hists_.find("diMuonLogger_")->second->getNbinsY()){
00539           // log runnumber, lumi block, event number & some
00540           // more pysics infomation for interesting events
00541           fill("diMuonLogger_", 0.5, diMuonLogged_+0.5, event.eventAuxiliary().run()); 
00542           fill("diMuonLogger_", 1.5, diMuonLogged_+0.5, event.eventAuxiliary().luminosityBlock()); 
00543           fill("diMuonLogger_", 2.5, diMuonLogged_+0.5, event.eventAuxiliary().event()); 
00544           fill("diMuonLogger_", 3.5, diMuonLogged_+0.5, isoMuons[0]->pt()); 
00545           fill("diMuonLogger_", 4.5, diMuonLogged_+0.5, isoMuons[1]->pt()); 
00546           if(leadingJets.size()>0) fill("diMuonLogger_", 5.5, diMuonLogged_+0.5, leadingJets[0].pt()); 
00547           if(leadingJets.size()>1) fill("diMuonLogger_", 6.5, diMuonLogged_+0.5, leadingJets[1].pt()); 
00548           fill("diMuonLogger_", 7.5, diMuonLogged_+0.5, caloMET.et()); 
00549           ++diMuonLogged_; 
00550         }
00551       }
00552     }
00553 
00554     // DIELEC channel
00555     if( decayChannel(isoMuons, isoElecs) == DIELEC ){
00556       fill("decayChannel_", 2.5);
00557       int charge = isoElecs[0]->charge()*isoElecs[1]->charge();
00558       double mass = (isoElecs[0]->p4()+isoElecs[1]->p4()).mass();
00559       fill(charge<0 ? "invMass_"    : "invMassWC_"    , mass       );
00560       fill(charge<0 ? "invMassLog_" : "invMassWCLog_" , log10(mass));
00561       if((lowerEdge_==-1. && upperEdge_==-1.) || (lowerEdge_<mass && mass<upperEdge_) ){
00562         fill("dEtaL1L2_"  , isoElecs[0]->eta()-isoElecs[1]->eta() );
00563         fill("sumEtaL1L2_", (isoElecs[0]->eta()+isoElecs[1]->eta())/2);
00564         fill("dPhiL1L2_"  , reco::deltaPhi(isoElecs[0]->phi(),isoElecs[1]->phi()) );
00565         fill("elecPt_", isoElecs[0]->pt()); fill("elecPt_", isoElecs[1]->pt()); 
00566         fill("lep1Pt_", isoElecs[0]->pt()); fill("lep2Pt_", isoElecs[1]->pt()); 
00567         if(diElecLogged_<=hists_.find("diElecLogger_")->second->getNbinsY()){
00568           // log runnumber, lumi block, event number & some
00569           // more pysics infomation for interesting events
00570           fill("diElecLogger_", 0.5, diElecLogged_+0.5, event.eventAuxiliary().run()); 
00571           fill("diElecLogger_", 1.5, diElecLogged_+0.5, event.eventAuxiliary().luminosityBlock()); 
00572           fill("diElecLogger_", 2.5, diElecLogged_+0.5, event.eventAuxiliary().event()); 
00573           fill("diElecLogger_", 3.5, diElecLogged_+0.5, isoElecs[0]->pt()); 
00574           fill("diElecLogger_", 4.5, diElecLogged_+0.5, isoElecs[1]->pt()); 
00575           if(leadingJets.size()>0) fill("diElecLogger_", 5.5, diElecLogged_+0.5, leadingJets[0].pt()); 
00576           if(leadingJets.size()>1) fill("diElecLogger_", 6.5, diElecLogged_+0.5, leadingJets[1].pt()); 
00577           fill("diElecLogger_", 7.5, diElecLogged_+0.5, caloMET.et()); 
00578           ++diElecLogged_; 
00579         }
00580       }
00581    }
00582   }
00583   
00584 }
00585 
00586 TopHLTDiLeptonOfflineDQM::TopHLTDiLeptonOfflineDQM(const edm::ParameterSet& cfg): triggerTable_(""), vertex_(""), vertexSelect_(0), beamspotSelect_(0)
00587 {
00588   // configure the preselection
00589   edm::ParameterSet presel=cfg.getParameter<edm::ParameterSet>("preselection");
00590   if( presel.existsAs<edm::ParameterSet>("trigger") ){
00591     edm::ParameterSet trigger=presel.getParameter<edm::ParameterSet>("trigger");
00592     triggerTable_=trigger.getParameter<edm::InputTag>("src");
00593     triggerPaths_=trigger.getParameter<std::vector<std::string> >("select");
00594   } 
00595   if( presel.existsAs<edm::ParameterSet>("vertex" ) ){
00596     edm::ParameterSet vertex=presel.getParameter<edm::ParameterSet>("vertex");
00597     vertex_= vertex.getParameter<edm::InputTag>("src");
00598     vertexSelect_= new StringCutObjectSelector<reco::Vertex>(vertex.getParameter<std::string>("select"));
00599   }
00600   if( presel.existsAs<edm::ParameterSet>("beamspot" ) ){
00601     edm::ParameterSet beamspot=presel.getParameter<edm::ParameterSet>("beamspot");
00602     beamspot_= beamspot.getParameter<edm::InputTag>("src");
00603     beamspotSelect_= new StringCutObjectSelector<reco::BeamSpot>(beamspot.getParameter<std::string>("select"));
00604   }
00605 
00606   // conifgure the selection
00607   std::vector<edm::ParameterSet> sel=cfg.getParameter<std::vector<edm::ParameterSet> >("selection");
00608   for(unsigned int i=0; i<sel.size(); ++i){
00609     selectionOrder_.push_back(sel.at(i).getParameter<std::string>("label"));
00610     selection_[selectionStep(selectionOrder_.back())] = std::make_pair(sel.at(i), new TopDiLeptonOffline::MonitorEnsemble(selectionStep(selectionOrder_.back()).c_str(), cfg.getParameter<edm::ParameterSet>("setup")));
00611   }
00612 }
00613 
00614 void 
00615 TopHLTDiLeptonOfflineDQM::analyze(const edm::Event& event, const edm::EventSetup& setup)
00616 { 
00617   if(!triggerTable_.label().empty()){
00618     edm::Handle<edm::TriggerResults> triggerTable;
00619     if( !event.getByLabel(triggerTable_, triggerTable) ) return;
00620     if(!acceptHLT(event, *triggerTable, triggerPaths_)) return;
00621   }
00622   if(!vertex_.label().empty()){
00623     edm::Handle<std::vector<reco::Vertex> > vertex;
00624     if( !event.getByLabel(vertex_, vertex) ) return;
00625     if(vertex->empty() || !(*vertexSelect_)(vertex->front())) return;
00626   }
00627   if(!beamspot_.label().empty()){
00628     edm::Handle<reco::BeamSpot> beamspot;
00629     if( !event.getByLabel(beamspot_, beamspot) ) return;
00630     if(!(*beamspotSelect_)(*beamspot)) return;
00631   }
00632   // apply selection steps
00633   for(std::vector<std::string>::const_iterator selIt=selectionOrder_.begin(); selIt!=selectionOrder_.end(); ++selIt){
00634     std::string key = selectionStep(*selIt), type = objectType(*selIt);
00635     if(selection_.find(key)!=selection_.end()){
00636       if(type=="empty"){
00637         selection_[key].second->fill(event, setup);
00638       }
00639       if(type=="Hlt" ){
00640 //      cout<<"HLT filled"<<endl;
00641         selection_[key].second->fill(event, setup);
00642       }
00643       if(type=="muons"){
00644         SelectionStepHLT<reco::Muon> step(selection_[key].first);
00645         if(step.select(event)){
00646           selection_[key].second->fill(event, setup);
00647         } else break;
00648       }
00649       if(type=="elecs"){
00650         SelectionStepHLT<reco::GsfElectron> step(selection_[key].first);
00651         if(step.select(event)){ 
00652           selection_[key].second->fill(event, setup);
00653         } else break;
00654       }
00655       if(type=="jets" ){
00656         SelectionStepHLT<reco::Jet> step(selection_[key].first);
00657         if(step.select(event, setup)){
00658           selection_[key].second->fill(event, setup);
00659         } else break;
00660       }
00661       if(type=="jets/pf" ){
00662         SelectionStepHLT<reco::PFJet> step(selection_[key].first);
00663         if(step.select(event, setup)){
00664           selection_[key].second->fill(event, setup);
00665         } else break;
00666       }
00667       if(type=="jets/calo" ){
00668         SelectionStepHLT<reco::CaloJet> step(selection_[key].first);
00669         if(step.select(event, setup)){
00670           selection_[key].second->fill(event, setup);
00671         } else break;
00672       }
00673       if(type=="met" ){
00674         SelectionStepHLT<reco::MET> step(selection_[key].first);
00675         if(step.select(event)){
00676           selection_[key].second->fill(event, setup);
00677         } else break;
00678       }
00679     }
00680   }
00681 }
00682 
00683