Reference Manual

Go to the documentation of this file.
00001 /*
00002  *  See header file for a description of this class.
00003  *
00004  *  $Date: 2010/03/25 11:03:38 $
00005  *  $Revision: 1.7 $
00006  *  \author A.Apresyan - Caltech
00007  */
00008 
00009 #include "DQMOffline/JetMET/interface/TcMETAnalyzer.h"
00010 #include "DataFormats/Common/interface/Handle.h"
00011 #include "FWCore/Common/interface/TriggerNames.h"
00012 
00013 #include "DataFormats/Math/interface/LorentzVector.h"
00014 
00015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00016 
00017 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
00018 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
00019 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
00020 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
00021 
00022 #include "DataFormats/Math/interface/LorentzVector.h"
00023 
00024 #include <string>
00025 using namespace edm;
00026 using namespace reco;
00027 using namespace math;
00028 
00029 // ***********************************************************
00030 TcMETAnalyzer::TcMETAnalyzer(const edm::ParameterSet& pSet) {
00031 
00032   parameters = pSet;
00033 
00034 }
00035 
00036 // ***********************************************************
00037 TcMETAnalyzer::~TcMETAnalyzer() { }
00038 
00039 void TcMETAnalyzer::beginJob(DQMStore * dbe) {
00040 
00041   evtCounter = 0;
00042   metname = "tcMETAnalyzer";
00043 
00044   // trigger information
00045   HLTPathsJetMBByName_ = parameters.getParameter<std::vector<std::string > >("HLTPathsJetMB");
00046 
00047   _hlt_HighPtJet = parameters.getParameter<std::string>("HLT_HighPtJet");
00048   _hlt_LowPtJet  = parameters.getParameter<std::string>("HLT_LowPtJet");
00049   _hlt_HighMET   = parameters.getParameter<std::string>("HLT_HighMET");
00050   _hlt_LowMET    = parameters.getParameter<std::string>("HLT_LowMET");
00051   _hlt_Ele       = parameters.getParameter<std::string>("HLT_Ele");
00052   _hlt_Muon      = parameters.getParameter<std::string>("HLT_Muon");
00053 
00054   // TcMET information
00055   theTcMETCollectionLabel       = parameters.getParameter<edm::InputTag>("TcMETCollectionLabel");
00056   _source                       = parameters.getParameter<std::string>("Source");
00057 
00058   // Other data collections
00059   HcalNoiseRBXCollectionTag   = parameters.getParameter<edm::InputTag>("HcalNoiseRBXCollection");
00060   HcalNoiseSummaryTag         = parameters.getParameter<edm::InputTag>("HcalNoiseSummary");
00061   theJetCollectionLabel       = parameters.getParameter<edm::InputTag>("JetCollectionLabel");
00062 
00063   // misc
00064   _verbose     = parameters.getParameter<int>("verbose");
00065   _etThreshold = parameters.getParameter<double>("etThreshold"); // MET threshold
00066   _allhist     = parameters.getParameter<bool>("allHist");       // Full set of monitoring histograms
00067   _allSelection= parameters.getParameter<bool>("allSelection");  // Plot with all sets of event selection
00068 
00069   _highPtTcJetThreshold = parameters.getParameter<double>("HighPtTcJetThreshold"); // High Pt Jet threshold
00070   _lowPtTcJetThreshold = parameters.getParameter<double>("LowPtTcJetThreshold");   // Low Pt Jet threshold
00071   _highTcMETThreshold = parameters.getParameter<double>("HighTcMETThreshold");     // High MET threshold
00072   _lowTcMETThreshold = parameters.getParameter<double>("LowTcMETThreshold");       // Low MET threshold
00073 
00074   //
00075   jetID = new reco::helper::JetIDHelper(parameters.getParameter<ParameterSet>("JetIDParams"));
00076 
00077   // DQStore stuff
00078   LogTrace(metname)<<"[TcMETAnalyzer] Parameters initialization";
00079   std::string DirName = "JetMET/MET/"+_source;
00080   dbe->setCurrentFolder(DirName);
00081 
00082   metME = dbe->book1D("metReco", "metReco", 4, 1, 5);
00083   metME->setBinLabel(2,"TcMET",1);
00084 
00085   _dbe = dbe;
00086 
00087   _FolderNames.push_back("All");
00088   _FolderNames.push_back("Cleanup");
00089   _FolderNames.push_back("HcalNoiseFilter");
00090   _FolderNames.push_back("HcalNoiseFilterTight");
00091   _FolderNames.push_back("JetID");
00092   _FolderNames.push_back("JetIDTight");
00093 
00094   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); 
00095        ic != _FolderNames.end(); ic++){
00096     if (*ic=="All")                  bookMESet(DirName+"/"+*ic);
00097     if (*ic=="Cleanup")              bookMESet(DirName+"/"+*ic);
00098     if (_allSelection){
00099     if (*ic=="HcalNoiseFilter")      bookMESet(DirName+"/"+*ic);
00100     if (*ic=="HcalNoiseFilterTight") bookMESet(DirName+"/"+*ic);
00101     if (*ic=="JetID")                bookMESet(DirName+"/"+*ic);
00102     if (*ic=="JetIDTight")           bookMESet(DirName+"/"+*ic);
00103     }
00104   }
00105 }
00106 
00107 // ***********************************************************
00108 void TcMETAnalyzer::endJob() {
00109 
00110   delete jetID;
00111 
00112 }
00113 
00114 // ***********************************************************
00115 void TcMETAnalyzer::bookMESet(std::string DirName)
00116 {
00117 
00118   bool bLumiSecPlot=false;
00119   if (DirName.find("All")!=std::string::npos) bLumiSecPlot=true;
00120 
00121   bookMonitorElement(DirName,bLumiSecPlot);
00122 
00123   if (_hlt_HighPtJet.size()){
00124     bookMonitorElement(DirName+"/"+"HighPtJet",false);
00125     meTriggerName_HighPtJet = _dbe->bookString("triggerName_HighPtJet", _hlt_HighPtJet);
00126   }  
00127 
00128   if (_hlt_LowPtJet.size()){
00129     bookMonitorElement(DirName+"/"+"LowPtJet",false);
00130     meTriggerName_LowPtJet = _dbe->bookString("triggerName_LowPtJet", _hlt_LowPtJet);
00131   }
00132 
00133   if (_hlt_HighMET.size()){
00134     bookMonitorElement(DirName+"/"+"HighMET",false);
00135     meTriggerName_HighMET = _dbe->bookString("triggerName_HighMET", _hlt_HighMET);
00136   }
00137 
00138   if (_hlt_LowMET.size()){
00139     bookMonitorElement(DirName+"/"+"LowMET",false);
00140     meTriggerName_LowMET = _dbe->bookString("triggerName_LowMET", _hlt_LowMET);
00141   }
00142 
00143   if (_hlt_Ele.size()){
00144     bookMonitorElement(DirName+"/"+"Ele",false);
00145     meTriggerName_Ele = _dbe->bookString("triggerName_Ele", _hlt_Ele);
00146   }
00147 
00148   if (_hlt_Muon.size()){
00149     bookMonitorElement(DirName+"/"+"Muon",false);
00150     meTriggerName_Muon = _dbe->bookString("triggerName_Muon", _hlt_Muon);
00151   }
00152 
00153 }
00154 
00155 // ***********************************************************
00156 void TcMETAnalyzer::bookMonitorElement(std::string DirName, bool bLumiSecPlot=false)
00157 {
00158 
00159   if (_verbose) std::cout << "booMonitorElement " << DirName << std::endl;
00160   _dbe->setCurrentFolder(DirName);
00161  
00162   meNevents              = _dbe->book1D("METTask_Nevents", "METTask_Nevents"   ,1,0,1);
00163   meTcMEx                = _dbe->book1D("METTask_TcMEx",   "METTask_TcMEx"   ,500,-500,500);
00164   meTcMEy                = _dbe->book1D("METTask_TcMEy",   "METTask_TcMEy"   ,500,-500,500);
00165   meTcEz                 = _dbe->book1D("METTask_TcEz",    "METTask_TcEz"    ,500,-500,500);
00166   meTcMETSig             = _dbe->book1D("METTask_TcMETSig","METTask_TcMETSig",51,0,51);
00167   meTcMET                = _dbe->book1D("METTask_TcMET",   "METTask_TcMET"   ,500,0,1000);
00168   meTcMETPhi             = _dbe->book1D("METTask_TcMETPhi","METTask_TcMETPhi",80,-TMath::Pi(),TMath::Pi());
00169   meTcSumET              = _dbe->book1D("METTask_TcSumET", "METTask_TcSumET" ,500,0,2000);
00170 
00171   meTcNeutralEMFraction  = _dbe->book1D("METTask_TcNeutralEMFraction", "METTask_TcNeutralEMFraction" ,50,0.,1.);
00172   meTcNeutralHadFraction = _dbe->book1D("METTask_TcNeutralHadFraction","METTask_TcNeutralHadFraction",50,0.,1.);
00173   meTcChargedEMFraction  = _dbe->book1D("METTask_TcChargedEMFraction", "METTask_TcChargedEMFraction" ,50,0.,1.);
00174   meTcChargedHadFraction = _dbe->book1D("METTask_TcChargedHadFraction","METTask_TcChargedHadFraction",50,0.,1.);
00175   meTcMuonFraction       = _dbe->book1D("METTask_TcMuonFraction",      "METTask_TcMuonFraction"      ,50,0.,1.);
00176 
00177   meTcMETIonFeedbck      = _dbe->book1D("METTask_TcMETIonFeedbck", "METTask_TcMETIonFeedbck" ,500,0,1000);
00178   meTcMETHPDNoise        = _dbe->book1D("METTask_TcMETHPDNoise",   "METTask_TcMETHPDNoise"   ,500,0,1000);
00179   meTcMETRBXNoise        = _dbe->book1D("METTask_TcMETRBXNoise",   "METTask_TcMETRBXNoise"   ,500,0,1000);
00180 
00181   if (_allhist){
00182     if (bLumiSecPlot){
00183       meTcMExLS              = _dbe->book2D("METTask_TcMEx_LS","METTask_TcMEx_LS",200,-200,200,50,0.,500.);
00184       meTcMEyLS              = _dbe->book2D("METTask_TcMEy_LS","METTask_TcMEy_LS",200,-200,200,50,0.,500.);
00185     }
00186   }
00187 }
00188 
00189 // ***********************************************************
00190 void TcMETAnalyzer::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup)
00191 {
00192 
00193 }
00194 
00195 // ***********************************************************
00196 void TcMETAnalyzer::endRun(const edm::Run& iRun, const edm::EventSetup& iSetup, DQMStore * dbe)
00197 {
00198   
00199   //
00200   //--- Check the time length of the Run from the lumi section plots
00201 
00202   std::string dirName = "JetMET/MET/"+_source+"/";
00203   _dbe->setCurrentFolder(dirName);
00204 
00205   TH1F* tlumisec;
00206 
00207   MonitorElement *meLumiSec = _dbe->get("aaa");
00208   meLumiSec = _dbe->get("JetMET/lumisec");
00209 
00210   int totlsec=0;
00211   double totltime=0.;
00212   if ( meLumiSec->getRootObject() ) {
00213     tlumisec = meLumiSec->getTH1F();
00214     for (int i=0; i<500; i++){
00215       if (tlumisec->GetBinContent(i+1)) totlsec++;
00216     }
00217     totltime = double(totlsec*90); // one lumi sec ~ 90 (sec)
00218   }
00219 
00220   if (totltime==0.) totltime=1.; 
00221 
00222   //
00223   //--- Make the integrated plots with rate (Hz)
00224 
00225   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); ic != _FolderNames.end(); ic++)
00226     {
00227 
00228       std::string DirName;
00229       DirName = dirName+*ic;
00230 
00231       makeRatePlot(DirName,totltime);
00232       if (_hlt_HighPtJet.size()) makeRatePlot(DirName+"/"+_hlt_HighPtJet,totltime);
00233       if (_hlt_LowPtJet.size())  makeRatePlot(DirName+"/"+_hlt_LowPtJet,totltime);
00234       if (_hlt_HighMET.size())   makeRatePlot(DirName+"/"+_hlt_HighMET,totltime);
00235       if (_hlt_LowMET.size())    makeRatePlot(DirName+"/"+_hlt_LowMET,totltime);
00236       if (_hlt_Ele.size())       makeRatePlot(DirName+"/"+_hlt_Ele,totltime);
00237       if (_hlt_Muon.size())      makeRatePlot(DirName+"/"+_hlt_Muon,totltime);
00238 
00239     }
00240 }
00241 
00242 // ***********************************************************
00243 void TcMETAnalyzer::makeRatePlot(std::string DirName, double totltime)
00244 {
00245 
00246   _dbe->setCurrentFolder(DirName);
00247   MonitorElement *meTcMET = _dbe->get(DirName+"/"+"METTask_TcMET");
00248 
00249   TH1F* tTcMET;
00250   TH1F* tTcMETRate;
00251 
00252   if ( meTcMET )
00253     if ( meTcMET->getRootObject() ) {
00254       tTcMET     = meTcMET->getTH1F();
00255       
00256       // Integral plot & convert number of events to rate (hz)
00257       tTcMETRate = (TH1F*) tTcMET->Clone("METTask_TcMETRate");
00258       for (int i = tTcMETRate->GetNbinsX()-1; i>=0; i--){
00259         tTcMETRate->SetBinContent(i+1,tTcMETRate->GetBinContent(i+2)+tTcMET->GetBinContent(i+1));
00260       }
00261       for (int i = 0; i<tTcMETRate->GetNbinsX(); i++){
00262         tTcMETRate->SetBinContent(i+1,tTcMETRate->GetBinContent(i+1)/double(totltime));
00263       }      
00264 
00265       meTcMETRate      = _dbe->book1D("METTask_TcMETRate",tTcMETRate);
00266       
00267     }
00268 
00269 }
00270 
00271 // ***********************************************************
00272 void TcMETAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup, 
00273                             const edm::TriggerResults& triggerResults) {
00274 
00275   if (_verbose) std::cout << "TcMETAnalyzer analyze" << std::endl;
00276 
00277   LogTrace(metname)<<"[TcMETAnalyzer] Analyze TcMET";
00278 
00279   metME->Fill(2);
00280 
00281   // ==========================================================  
00282   // Trigger information 
00283   //
00284   _trig_JetMB=0;
00285   _trig_HighPtJet=0;
00286   _trig_LowPtJet=0;
00287   _trig_HighMET=0;
00288   _trig_LowMET=0;
00289   if(&triggerResults) {   
00290 
00292 
00293     //
00294     //
00295     // Check how many HLT triggers are in triggerResults 
00296     int ntrigs = triggerResults.size();
00297     if (_verbose) std::cout << "ntrigs=" << ntrigs << std::endl;
00298     
00299     //
00300     //
00301     // If index=ntrigs, this HLT trigger doesn't exist in the HLT table for this data.
00302     const edm::TriggerNames & triggerNames = iEvent.triggerNames(triggerResults);
00303     
00304     //
00305     //
00306     // count number of requested Jet or MB HLT paths which have fired
00307     for (unsigned int i=0; i!=HLTPathsJetMBByName_.size(); i++) {
00308       unsigned int triggerIndex = triggerNames.triggerIndex(HLTPathsJetMBByName_[i]);
00309       if (triggerIndex<triggerResults.size()) {
00310         if (triggerResults.accept(triggerIndex)) {
00311           _trig_JetMB++;
00312         }
00313       }
00314     }
00315     // for empty input vectors (n==0), take all HLT triggers!
00316     if (HLTPathsJetMBByName_.size()==0) _trig_JetMB=triggerResults.size()-1;
00317 
00318     //
00319     if (_verbose) std::cout << "triggerNames size" << " " << triggerNames.size() << std::endl;
00320     if (_verbose) std::cout << _hlt_HighPtJet << " " << triggerNames.triggerIndex(_hlt_HighPtJet) << std::endl;
00321     if (_verbose) std::cout << _hlt_LowPtJet  << " " << triggerNames.triggerIndex(_hlt_LowPtJet)  << std::endl;
00322     if (_verbose) std::cout << _hlt_HighMET   << " " << triggerNames.triggerIndex(_hlt_HighMET)   << std::endl;
00323     if (_verbose) std::cout << _hlt_LowMET    << " " << triggerNames.triggerIndex(_hlt_LowMET)    << std::endl;
00324     if (_verbose) std::cout << _hlt_Ele       << " " << triggerNames.triggerIndex(_hlt_Ele)       << std::endl;
00325     if (_verbose) std::cout << _hlt_Muon      << " " << triggerNames.triggerIndex(_hlt_Muon)      << std::endl;
00326 
00327     if (triggerNames.triggerIndex(_hlt_HighPtJet) != triggerNames.size() &&
00328         triggerResults.accept(triggerNames.triggerIndex(_hlt_HighPtJet))) _trig_HighPtJet=1;
00329 
00330     if (triggerNames.triggerIndex(_hlt_LowPtJet)  != triggerNames.size() &&
00331         triggerResults.accept(triggerNames.triggerIndex(_hlt_LowPtJet)))  _trig_LowPtJet=1;
00332 
00333     if (triggerNames.triggerIndex(_hlt_HighMET)   != triggerNames.size() &&
00334         triggerResults.accept(triggerNames.triggerIndex(_hlt_HighMET)))   _trig_HighMET=1;
00335 
00336     if (triggerNames.triggerIndex(_hlt_LowMET)    != triggerNames.size() &&
00337         triggerResults.accept(triggerNames.triggerIndex(_hlt_LowMET)))    _trig_LowMET=1;
00338 
00339     if (triggerNames.triggerIndex(_hlt_Ele)       != triggerNames.size() &&
00340         triggerResults.accept(triggerNames.triggerIndex(_hlt_Ele)))       _trig_Ele=1;
00341 
00342     if (triggerNames.triggerIndex(_hlt_Muon)      != triggerNames.size() &&
00343         triggerResults.accept(triggerNames.triggerIndex(_hlt_Muon)))      _trig_Muon=1;
00344     
00345   } else {
00346 
00347     edm::LogInfo("TcMetAnalyzer") << "TriggerResults::HLT not found, "
00348       "automatically select events"; 
00349 
00350     // TriggerResults object not found. Look at all events.    
00351     _trig_JetMB=1;
00352   }
00353 
00354   // ==========================================================
00355   // TcMET information
00356   
00357   // **** Get the MET container  
00358   edm::Handle<reco::METCollection> tcmetcoll;
00359   iEvent.getByLabel(theTcMETCollectionLabel, tcmetcoll);
00360   
00361   if(!tcmetcoll.isValid()) return;
00362 
00363   const METCollection *tcmetcol = tcmetcoll.product();
00364   const MET *tcmet;
00365   tcmet = &(tcmetcol->front());
00366     
00367   LogTrace(metname)<<"[TcMETAnalyzer] Call to the TcMET analyzer";
00368 
00369   // ==========================================================
00370   //
00371   edm::Handle<HcalNoiseRBXCollection> HRBXCollection;
00372   iEvent.getByLabel(HcalNoiseRBXCollectionTag,HRBXCollection);
00373   if (!HRBXCollection.isValid()) {
00374     LogDebug("") << "TcMETAnalyzer: Could not find HcalNoiseRBX Collection" << std::endl;
00375     if (_verbose) std::cout << "TcMETAnalyzer: Could not find HcalNoiseRBX Collection" << std::endl;
00376   }
00377   
00378   edm::Handle<HcalNoiseSummary> HNoiseSummary;
00379   iEvent.getByLabel(HcalNoiseSummaryTag,HNoiseSummary);
00380   if (!HNoiseSummary.isValid()) {
00381     LogDebug("") << "TcMETAnalyzer: Could not find Hcal NoiseSummary product" << std::endl;
00382     if (_verbose) std::cout << "TcMETAnalyzer: Could not find Hcal NoiseSummary product" << std::endl;
00383   }
00384 
00385   edm::Handle<reco::CaloJetCollection> caloJets;
00386   iEvent.getByLabel(theJetCollectionLabel, caloJets);
00387   if (!caloJets.isValid()) {
00388     LogDebug("") << "TcMETAnalyzer: Could not find jet product" << std::endl;
00389     if (_verbose) std::cout << "TcMETAnalyzer: Could not find jet product" << std::endl;
00390   }
00391 
00392   // ==========================================================
00393   // TcMET sanity check
00394 
00395   //   if (_source=="TcMET") validateMET(*tcmet, tcCandidates);
00396   
00397   // ==========================================================
00398   // JetID 
00399 
00400   if (_verbose) std::cout << "JetID starts" << std::endl;
00401   
00402   //
00403   // --- Loose cuts, not Tc specific for now!
00404   //
00405   bool bJetID=true;
00406   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00407        cal!=caloJets->end(); ++cal){ 
00408     jetID->calculate(iEvent, *cal);
00409     if (_verbose) std::cout << jetID->n90Hits() << " " 
00410                             << jetID->restrictedEMF() << " "
00411                             << cal->pt() << std::endl;
00412     if (cal->pt()>10.){
00413       //
00414       // for all regions
00415       if (jetID->n90Hits()<2)  bJetID=false; 
00416       if (jetID->fHPD()>=0.98) bJetID=false; 
00417       //if (jetID->restrictedEMF()<0.01) bJetID=false; 
00418       //
00419       // for non-forward
00420       if (fabs(cal->eta())<2.55){
00421         if (cal->emEnergyFraction()<=0.01) bJetID=false; 
00422       }
00423       // for forward
00424       else {
00425         if (cal->emEnergyFraction()<=-0.9) bJetID=false; 
00426         if (cal->pt()>80.){
00427         if (cal->emEnergyFraction()>= 1.0) bJetID=false; 
00428         }
00429       } // forward vs non-forward
00430     }   // pt>10 GeV/c
00431   }     // calor-jets loop
00432  
00433   //
00434   // --- Tight cuts
00435   //
00436   bool bJetIDTight=true;
00437   bJetIDTight=bJetID;
00438   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00439        cal!=caloJets->end(); ++cal){
00440     jetID->calculate(iEvent, *cal);
00441     if (cal->pt()>25.){
00442       //
00443       // for all regions
00444       if (jetID->fHPD()>=0.95) bJetIDTight=false; 
00445       //
00446       // for 1.0<|eta|<1.75
00447       if (fabs(cal->eta())>=1.00 && fabs(cal->eta())<1.75){
00448         if (cal->pt()>80. && cal->emEnergyFraction()>=1.) bJetIDTight=false; 
00449       }
00450       //
00451       // for 1.75<|eta|<2.55
00452       else if (fabs(cal->eta())>=1.75 && fabs(cal->eta())<2.55){
00453         if (cal->pt()>80. && cal->emEnergyFraction()>=1.) bJetIDTight=false; 
00454       }
00455       //
00456       // for 2.55<|eta|<3.25
00457       else if (fabs(cal->eta())>=2.55 && fabs(cal->eta())<3.25){
00458         if (cal->pt()< 50.                   && cal->emEnergyFraction()<=-0.3) bJetIDTight=false; 
00459         if (cal->pt()>=50. && cal->pt()< 80. && cal->emEnergyFraction()<=-0.2) bJetIDTight=false; 
00460         if (cal->pt()>=80. && cal->pt()<340. && cal->emEnergyFraction()<=-0.1) bJetIDTight=false; 
00461         if (cal->pt()>=340.                  && cal->emEnergyFraction()<=-0.1 
00462                                              && cal->emEnergyFraction()>=0.95) bJetIDTight=false; 
00463       }
00464       //
00465       // for 3.25<|eta|
00466       else if (fabs(cal->eta())>=3.25){
00467         if (cal->pt()< 50.                   && cal->emEnergyFraction()<=-0.3
00468                                              && cal->emEnergyFraction()>=0.90) bJetIDTight=false; 
00469         if (cal->pt()>=50. && cal->pt()<130. && cal->emEnergyFraction()<=-0.2
00470                                              && cal->emEnergyFraction()>=0.80) bJetIDTight=false; 
00471         if (cal->pt()>=130.                  && cal->emEnergyFraction()<=-0.1 
00472                                              && cal->emEnergyFraction()>=0.70) bJetIDTight=false; 
00473       }
00474     }   // pt>10 GeV/c
00475   }     // calor-jets loop
00476   
00477   if (_verbose) std::cout << "JetID ends" << std::endl;
00478 
00479 
00480   // ==========================================================
00481   // HCAL Noise filter
00482   
00483   bool bHcalNoiseFilter      = HNoiseSummary->passLooseNoiseFilter();
00484   bool bHcalNoiseFilterTight = HNoiseSummary->passTightNoiseFilter();
00485 
00486   // ==========================================================
00487   // Reconstructed MET Information - fill MonitorElements
00488   
00489   std::string DirName = "JetMET/MET/"+_source;
00490   
00491   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); 
00492        ic != _FolderNames.end(); ic++){
00493     if (*ic=="All")                                   fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00494     if (*ic=="Cleanup" && bHcalNoiseFilter && bJetID) fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00495     if (_allSelection) {
00496     if (*ic=="HcalNoiseFilter"      && bHcalNoiseFilter )       fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00497     if (*ic=="HcalNoiseFilterTight" && bHcalNoiseFilterTight )  fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00498     if (*ic=="JetID"      && bJetID)                            fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00499     if (*ic=="JetIDTight" && bJetIDTight)                       fillMESet(iEvent, DirName+"/"+*ic, *tcmet);
00500     }
00501   }
00502 }
00503 
00504 // ***********************************************************
00505 void TcMETAnalyzer::fillMESet(const edm::Event& iEvent, std::string DirName, 
00506                               const reco::MET& tcmet)
00507 {
00508 
00509   _dbe->setCurrentFolder(DirName);
00510 
00511   bool bLumiSecPlot=false;
00512   if (DirName.find("All")) bLumiSecPlot=true;
00513 
00514   if (_trig_JetMB) fillMonitorElement(iEvent,DirName,"",tcmet, bLumiSecPlot);
00515   if (_hlt_HighPtJet.size() && _trig_HighPtJet) fillMonitorElement(iEvent,DirName,"HighPtJet",tcmet,false);
00516   if (_hlt_LowPtJet.size() && _trig_LowPtJet) fillMonitorElement(iEvent,DirName,"LowPtJet",tcmet,false);
00517   if (_hlt_HighMET.size() && _trig_HighMET) fillMonitorElement(iEvent,DirName,"HighMET",tcmet,false);
00518   if (_hlt_LowMET.size() && _trig_LowMET) fillMonitorElement(iEvent,DirName,"LowMET",tcmet,false);
00519   if (_hlt_Ele.size() && _trig_Ele) fillMonitorElement(iEvent,DirName,"Ele",tcmet,false);
00520   if (_hlt_Muon.size() && _trig_Muon) fillMonitorElement(iEvent,DirName,"Muon",tcmet,false);
00521 }
00522 
00523 // ***********************************************************
00524 void TcMETAnalyzer::fillMonitorElement(const edm::Event& iEvent, std::string DirName, 
00525                                          std::string TriggerTypeName, 
00526                                          const reco::MET& tcmet, bool bLumiSecPlot)
00527 {
00528 
00529   if (TriggerTypeName=="HighPtJet") {
00530     if (!selectHighPtJetEvent(iEvent)) return;
00531   }
00532   else if (TriggerTypeName=="LowPtJet") {
00533     if (!selectLowPtJetEvent(iEvent)) return;
00534   }
00535   else if (TriggerTypeName=="HighMET") {
00536     if (tcmet.pt()<_highTcMETThreshold) return;
00537   }
00538   else if (TriggerTypeName=="LowMET") {
00539     if (tcmet.pt()<_lowTcMETThreshold) return;
00540   }
00541   else if (TriggerTypeName=="Ele") {
00542     if (!selectWElectronEvent(iEvent)) return;
00543   }
00544   else if (TriggerTypeName=="Muon") {
00545     if (!selectWMuonEvent(iEvent)) return;
00546   }
00547   
00548 // Reconstructed MET Information
00549   double tcSumET  = tcmet.sumEt();
00550   double tcMETSig = tcmet.mEtSig();
00551   double tcEz     = tcmet.e_longitudinal();
00552   double tcMET    = tcmet.pt();
00553   double tcMEx    = tcmet.px();
00554   double tcMEy    = tcmet.py();
00555   double tcMETPhi = tcmet.phi();
00556 
00557   //
00558   int myLuminosityBlock;
00559   //  myLuminosityBlock = (evtCounter++)/1000;
00560   myLuminosityBlock = iEvent.luminosityBlock();
00561   //
00562 
00563   if (TriggerTypeName!="") DirName = DirName +"/"+TriggerTypeName;
00564 
00565   if (_verbose) std::cout << "_etThreshold = " << _etThreshold << std::endl;
00566   if (tcMET>_etThreshold){
00567     
00568     meTcMEx    = _dbe->get(DirName+"/"+"METTask_TcMEx");    if (meTcMEx    && meTcMEx->getRootObject())    meTcMEx->Fill(tcMEx);
00569     meTcMEy    = _dbe->get(DirName+"/"+"METTask_TcMEy");    if (meTcMEy    && meTcMEy->getRootObject())    meTcMEy->Fill(tcMEy);
00570     meTcMET    = _dbe->get(DirName+"/"+"METTask_TcMET");    if (meTcMET    && meTcMET->getRootObject())    meTcMET->Fill(tcMET);
00571     meTcMETPhi = _dbe->get(DirName+"/"+"METTask_TcMETPhi"); if (meTcMETPhi && meTcMETPhi->getRootObject()) meTcMETPhi->Fill(tcMETPhi);
00572     meTcSumET  = _dbe->get(DirName+"/"+"METTask_TcSumET");  if (meTcSumET  && meTcSumET->getRootObject())  meTcSumET->Fill(tcSumET);
00573     meTcMETSig = _dbe->get(DirName+"/"+"METTask_TcMETSig"); if (meTcMETSig && meTcMETSig->getRootObject()) meTcMETSig->Fill(tcMETSig);
00574     meTcEz     = _dbe->get(DirName+"/"+"METTask_TcEz");     if (meTcEz     && meTcEz->getRootObject())     meTcEz->Fill(tcEz);
00575 
00576     meTcMETIonFeedbck = _dbe->get(DirName+"/"+"METTask_TcMETIonFeedbck");  if (meTcMETIonFeedbck && meTcMETIonFeedbck->getRootObject()) meTcMETIonFeedbck->Fill(tcMET);
00577     meTcMETHPDNoise   = _dbe->get(DirName+"/"+"METTask_TcMETHPDNoise");    if (meTcMETHPDNoise   && meTcMETHPDNoise->getRootObject())   meTcMETHPDNoise->Fill(tcMET);
00578     meTcMETRBXNoise   = _dbe->get(DirName+"/"+"METTask_TcMETRBXNoise");    if (meTcMETRBXNoise   && meTcMETRBXNoise->getRootObject())   meTcMETRBXNoise->Fill(tcMET);
00579         
00580     if (_allhist){
00581       if (bLumiSecPlot){
00582         meTcMExLS = _dbe->get(DirName+"/"+"METTask_TcMExLS"); if (meTcMExLS && meTcMExLS->getRootObject()) meTcMExLS->Fill(tcMEx,myLuminosityBlock);
00583         meTcMEyLS = _dbe->get(DirName+"/"+"METTask_TcMEyLS"); if (meTcMEyLS && meTcMEyLS->getRootObject()) meTcMEyLS->Fill(tcMEy,myLuminosityBlock);
00584       }
00585     } // _allhist
00586   } // et threshold cut
00587 }
00588 
00589 // ***********************************************************
00590 bool TcMETAnalyzer::selectHighPtJetEvent(const edm::Event& iEvent){
00591 
00592   bool return_value=false;
00593 
00594   edm::Handle<reco::CaloJetCollection> caloJets;
00595   iEvent.getByLabel(theJetCollectionLabel, caloJets);
00596   if (!caloJets.isValid()) {
00597     LogDebug("") << "TcMETAnalyzer: Could not find jet product" << std::endl;
00598     if (_verbose) std::cout << "TcMETAnalyzer: Could not find jet product" << std::endl;
00599   }
00600 
00601   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00602        cal!=caloJets->end(); ++cal){
00603     if (cal->pt()>_highPtTcJetThreshold){
00604       return_value=true;
00605     }
00606   }
00607   
00608   return return_value;
00609 }
00610 
00611 // // ***********************************************************
00612 bool TcMETAnalyzer::selectLowPtJetEvent(const edm::Event& iEvent){
00613 
00614   bool return_value=false;
00615 
00616   edm::Handle<reco::CaloJetCollection> caloJets;
00617   iEvent.getByLabel(theJetCollectionLabel, caloJets);
00618   if (!caloJets.isValid()) {
00619     LogDebug("") << "TcMETAnalyzer: Could not find jet product" << std::endl;
00620     if (_verbose) std::cout << "TcMETAnalyzer: Could not find jet product" << std::endl;
00621   }
00622 
00623   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00624        cal!=caloJets->end(); ++cal){
00625     if (cal->pt()>_lowPtTcJetThreshold){
00626       return_value=true;
00627     }
00628   }
00629 
00630   return return_value;
00631 
00632 }
00633 
00634 // ***********************************************************
00635 bool TcMETAnalyzer::selectWElectronEvent(const edm::Event& iEvent){
00636 
00637   bool return_value=false;
00638 
00639   /*
00640     W-electron event selection comes here
00641    */
00642 
00643   return return_value;
00644 
00645 }
00646 
00647 // ***********************************************************
00648 bool TcMETAnalyzer::selectWMuonEvent(const edm::Event& iEvent){
00649 
00650   bool return_value=false;
00651 
00652   /*
00653     W-muon event selection comes here
00654    */
00655 
00656   return return_value;
00657 
00658 }
/data/refman/pasoursint/CMSSW_5_3_0/src/DQMOffline/JetMET/src/TcMETAnalyzer.cc
