Reference Manual

Go to the documentation of this file.
00001 /*
00002  *  See header file for a description of this class.
00003  *
00004  *  $Date: 2011/10/10 14:43:31 $
00005  *  $Revision: 1.44 $
00006  *  \author A.Apresyan - Caltech
00007  *          K.Hatakeyama - Baylor
00008  */
00009 
00010 #include "DQMOffline/JetMET/interface/METAnalyzer.h"
00011 #include "DataFormats/Common/interface/Handle.h"
00012 #include "FWCore/Common/interface/TriggerNames.h"
00013 
00014 #include "DataFormats/Math/interface/LorentzVector.h"
00015 
00016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00017 
00018 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
00019 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
00020 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
00021 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
00022 
00023 #include "DataFormats/Math/interface/LorentzVector.h"
00024 
00025 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h"
00026 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetup.h"
00027 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
00028 
00029 #include <string>
00030 using namespace edm;
00031 using namespace reco;
00032 using namespace math;
00033 
00034 // ***********************************************************
00035 METAnalyzer::METAnalyzer(const edm::ParameterSet& pSet) {
00036 
00037   parameters = pSet;
00038 
00039   edm::ParameterSet highptjetparms = parameters.getParameter<edm::ParameterSet>("highPtJetTrigger");
00040   edm::ParameterSet lowptjetparms  = parameters.getParameter<edm::ParameterSet>("lowPtJetTrigger" );
00041   edm::ParameterSet minbiasparms   = parameters.getParameter<edm::ParameterSet>("minBiasTrigger"  );
00042   edm::ParameterSet highmetparms   = parameters.getParameter<edm::ParameterSet>("highMETTrigger"  );
00043   edm::ParameterSet lowmetparms    = parameters.getParameter<edm::ParameterSet>("lowMETTrigger"   );
00044   edm::ParameterSet eleparms       = parameters.getParameter<edm::ParameterSet>("eleTrigger"      );
00045   edm::ParameterSet muonparms      = parameters.getParameter<edm::ParameterSet>("muonTrigger"     );
00046 
00047   //genericTriggerEventFlag_( new GenericTriggerEventFlag( conf_ ) );
00048   _HighPtJetEventFlag = new GenericTriggerEventFlag( highptjetparms );
00049   _LowPtJetEventFlag  = new GenericTriggerEventFlag( lowptjetparms  );
00050   _MinBiasEventFlag   = new GenericTriggerEventFlag( minbiasparms   );
00051   _HighMETEventFlag   = new GenericTriggerEventFlag( highmetparms   );
00052   _LowMETEventFlag    = new GenericTriggerEventFlag( lowmetparms    );
00053   _EleEventFlag       = new GenericTriggerEventFlag( eleparms       );
00054   _MuonEventFlag      = new GenericTriggerEventFlag( muonparms      );
00055 
00056   highPtJetExpr_ = highptjetparms.getParameter<std::vector<std::string> >("hltPaths");
00057   lowPtJetExpr_  = lowptjetparms .getParameter<std::vector<std::string> >("hltPaths");
00058   highMETExpr_   = highmetparms  .getParameter<std::vector<std::string> >("hltPaths");
00059   lowMETExpr_    = lowmetparms   .getParameter<std::vector<std::string> >("hltPaths");
00060   muonExpr_      = muonparms     .getParameter<std::vector<std::string> >("hltPaths");
00061   elecExpr_      = eleparms      .getParameter<std::vector<std::string> >("hltPaths");
00062   minbiasExpr_   = minbiasparms  .getParameter<std::vector<std::string> >("hltPaths");
00063 
00064 }
00065 
00066 // ***********************************************************
00067 METAnalyzer::~METAnalyzer() { 
00068 
00069   delete _HighPtJetEventFlag;
00070   delete _LowPtJetEventFlag;
00071   delete _MinBiasEventFlag;
00072   delete _HighMETEventFlag;
00073   delete _LowMETEventFlag;
00074   delete _EleEventFlag;
00075   delete _MuonEventFlag;
00076 
00077 }
00078 
00079 void METAnalyzer::beginJob(DQMStore * dbe) {
00080 
00081   evtCounter = 0;
00082   metname = "METAnalyzer";
00083 
00084   // trigger information
00085   HLTPathsJetMBByName_ = parameters.getParameter<std::vector<std::string > >("HLTPathsJetMB");
00086 
00087   theCleaningParameters = parameters.getParameter<ParameterSet>("CleaningParameters"),
00088 
00089   //Trigger parameters
00090   gtTag          = theCleaningParameters.getParameter<edm::InputTag>("gtLabel");
00091   _techTrigsAND  = theCleaningParameters.getParameter<std::vector<unsigned > >("techTrigsAND");
00092   _techTrigsOR   = theCleaningParameters.getParameter<std::vector<unsigned > >("techTrigsOR");
00093   _techTrigsNOT  = theCleaningParameters.getParameter<std::vector<unsigned > >("techTrigsNOT");
00094 
00095   _doHLTPhysicsOn = theCleaningParameters.getParameter<bool>("doHLTPhysicsOn");
00096   _hlt_PhysDec    = theCleaningParameters.getParameter<std::string>("HLT_PhysDec");
00097 
00098   _tightBHFiltering     = theCleaningParameters.getParameter<bool>("tightBHFiltering");
00099   _tightJetIDFiltering  = theCleaningParameters.getParameter<int>("tightJetIDFiltering");
00100   _tightHcalFiltering   = theCleaningParameters.getParameter<bool>("tightHcalFiltering");
00101 
00102   // ==========================================================
00103   //DCS information
00104   // ==========================================================
00105   DCSFilter = new JetMETDQMDCSFilter(parameters.getParameter<ParameterSet>("DCSFilter"));
00106 
00107   //Vertex requirements
00108   _doPVCheck          = theCleaningParameters.getParameter<bool>("doPrimaryVertexCheck");
00109   vertexTag  = theCleaningParameters.getParameter<edm::InputTag>("vertexLabel");
00110 
00111   if (_doPVCheck) {
00112     _nvtx_min        = theCleaningParameters.getParameter<int>("nvtx_min");
00113     _nvtxtrks_min    = theCleaningParameters.getParameter<int>("nvtxtrks_min");
00114     _vtxndof_min     = theCleaningParameters.getParameter<int>("vtxndof_min");
00115     _vtxchi2_max     = theCleaningParameters.getParameter<double>("vtxchi2_max");
00116     _vtxz_max        = theCleaningParameters.getParameter<double>("vtxz_max");
00117   }
00118 
00119   // MET information
00120   theMETCollectionLabel       = parameters.getParameter<edm::InputTag>("METCollectionLabel");
00121   _source                     = parameters.getParameter<std::string>("Source");
00122 
00123   if (theMETCollectionLabel.label() == "tcMet" ) {
00124     inputTrackLabel         = parameters.getParameter<edm::InputTag>("InputTrackLabel");    
00125     inputMuonLabel          = parameters.getParameter<edm::InputTag>("InputMuonLabel");
00126     inputElectronLabel      = parameters.getParameter<edm::InputTag>("InputElectronLabel");
00127     inputBeamSpotLabel      = parameters.getParameter<edm::InputTag>("InputBeamSpotLabel");
00128   }
00129 
00130   // Other data collections
00131   theJetCollectionLabel       = parameters.getParameter<edm::InputTag>("JetCollectionLabel");
00132   HcalNoiseRBXCollectionTag   = parameters.getParameter<edm::InputTag>("HcalNoiseRBXCollection");
00133   HcalNoiseSummaryTag         = parameters.getParameter<edm::InputTag>("HcalNoiseSummary");
00134   BeamHaloSummaryTag          = parameters.getParameter<edm::InputTag>("BeamHaloSummaryLabel");
00135 
00136   // misc
00137   _verbose      = parameters.getParameter<int>("verbose");
00138   _etThreshold  = parameters.getParameter<double>("etThreshold"); // MET threshold
00139   _allhist      = parameters.getParameter<bool>("allHist");       // Full set of monitoring histograms
00140   _allSelection = parameters.getParameter<bool>("allSelection");  // Plot with all sets of event selection
00141   _cleanupSelection = parameters.getParameter<bool>("cleanupSelection");  // Plot with all sets of event selection
00142 
00143   _FolderName              = parameters.getUntrackedParameter<std::string>("FolderName");
00144 
00145   _highPtJetThreshold = parameters.getParameter<double>("HighPtJetThreshold"); // High Pt Jet threshold
00146   _lowPtJetThreshold  = parameters.getParameter<double>("LowPtJetThreshold");   // Low Pt Jet threshold
00147   _highMETThreshold   = parameters.getParameter<double>("HighMETThreshold");     // High MET threshold
00148   _lowMETThreshold    = parameters.getParameter<double>("LowMETThreshold");       // Low MET threshold
00149 
00150   //
00151   jetID = new reco::helper::JetIDHelper(parameters.getParameter<ParameterSet>("JetIDParams"));
00152 
00153   // DQStore stuff
00154   LogTrace(metname)<<"[METAnalyzer] Parameters initialization";
00155   std::string DirName = _FolderName+_source;
00156   dbe->setCurrentFolder(DirName);
00157 
00158   hmetME = dbe->book1D("metReco", "metReco", 4, 1, 5);
00159   hmetME->setBinLabel(2,"MET",1);
00160 
00161   _dbe = dbe;
00162 
00163   _FolderNames.push_back("All");
00164   _FolderNames.push_back("BasicCleanup");
00165   _FolderNames.push_back("ExtraCleanup");
00166   _FolderNames.push_back("HcalNoiseFilter");
00167   _FolderNames.push_back("HcalNoiseFilterTight");
00168   _FolderNames.push_back("JetIDMinimal");
00169   _FolderNames.push_back("JetIDLoose");
00170   _FolderNames.push_back("JetIDTight");
00171   _FolderNames.push_back("BeamHaloIDTightPass");
00172   _FolderNames.push_back("BeamHaloIDLoosePass");
00173   _FolderNames.push_back("Triggers");
00174   _FolderNames.push_back("PV");
00175 
00176   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); 
00177        ic != _FolderNames.end(); ic++){
00178     if (*ic=="All")                  bookMESet(DirName+"/"+*ic);
00179     if (_cleanupSelection){
00180     if (*ic=="BasicCleanup")         bookMESet(DirName+"/"+*ic);
00181     if (*ic=="ExtraCleanup")         bookMESet(DirName+"/"+*ic);
00182     }
00183     if (_allSelection){
00184       if (*ic=="HcalNoiseFilter")      bookMESet(DirName+"/"+*ic);
00185       if (*ic=="HcalNoiseFilterTight") bookMESet(DirName+"/"+*ic);
00186       if (*ic=="JetIDMinimal")         bookMESet(DirName+"/"+*ic);
00187       if (*ic=="JetIDLoose")           bookMESet(DirName+"/"+*ic);
00188       if (*ic=="JetIDTight")           bookMESet(DirName+"/"+*ic);
00189       if (*ic=="BeamHaloIDTightPass")  bookMESet(DirName+"/"+*ic);
00190       if (*ic=="BeamHaloIDLoosePass")  bookMESet(DirName+"/"+*ic);
00191       if (*ic=="Triggers")             bookMESet(DirName+"/"+*ic);
00192       if (*ic=="PV")                   bookMESet(DirName+"/"+*ic);
00193     }
00194   }
00195 }
00196 
00197 // ***********************************************************
00198 void METAnalyzer::endJob() {
00199 
00200   delete jetID;
00201   delete DCSFilter;
00202 
00203 }
00204 
00205 // ***********************************************************
00206 void METAnalyzer::bookMESet(std::string DirName)
00207 {
00208 
00209   bool bLumiSecPlot=false;
00210   if (DirName.find("All")!=std::string::npos) bLumiSecPlot=true;
00211 
00212   bookMonitorElement(DirName,bLumiSecPlot);
00213 
00214   if ( _HighPtJetEventFlag->on() ) {
00215     bookMonitorElement(DirName+"/"+"HighPtJet",false);
00216     hTriggerName_HighPtJet = _dbe->bookString("triggerName_HighPtJet", highPtJetExpr_[0]);
00217   }  
00218 
00219   if ( _LowPtJetEventFlag->on() ) {
00220     bookMonitorElement(DirName+"/"+"LowPtJet",false);
00221     hTriggerName_LowPtJet = _dbe->bookString("triggerName_LowPtJet", lowPtJetExpr_[0]);
00222   }
00223 
00224   if ( _MinBiasEventFlag->on() ) {
00225     bookMonitorElement(DirName+"/"+"MinBias",false);
00226     hTriggerName_MinBias = _dbe->bookString("triggerName_MinBias", minbiasExpr_[0]);
00227     if (_verbose) std::cout << "_MinBiasEventFlag is on, folder created\n";
00228   }
00229 
00230   if ( _HighMETEventFlag->on() ) {
00231     bookMonitorElement(DirName+"/"+"HighMET",false);
00232     hTriggerName_HighMET = _dbe->bookString("triggerName_HighMET", highMETExpr_[0]);
00233   }
00234 
00235   if ( _LowMETEventFlag->on() ) {
00236     bookMonitorElement(DirName+"/"+"LowMET",false);
00237     hTriggerName_LowMET = _dbe->bookString("triggerName_LowMET", lowMETExpr_[0]);
00238   }
00239 
00240   if ( _EleEventFlag->on() ) {
00241     bookMonitorElement(DirName+"/"+"Ele",false);
00242     hTriggerName_Ele = _dbe->bookString("triggerName_Ele", elecExpr_[0]);
00243     if (_verbose) std::cout << "_EleEventFlag is on, folder created\n";
00244   }
00245 
00246   if ( _MuonEventFlag->on() ) {
00247     bookMonitorElement(DirName+"/"+"Muon",false);
00248     hTriggerName_Muon = _dbe->bookString("triggerName_Muon", muonExpr_[0]);
00249     if (_verbose) std::cout << "_MuonEventFlag is on, folder created\n";
00250   }
00251 }
00252 
00253 // ***********************************************************
00254 void METAnalyzer::bookMonitorElement(std::string DirName, bool bLumiSecPlot=false)
00255 {
00256 
00257   if (_verbose) std::cout << "bookMonitorElement " << DirName << std::endl;
00258   _dbe->setCurrentFolder(DirName);
00259  
00260   //hNevents            = _dbe->book1D("METTask_Nevents", "METTask_Nevents"   ,1,0,1);
00261   hMEx                = _dbe->book1D("METTask_MEx",   "METTask_MEx"   ,200,-500,500);
00262   hMEx->setAxisTitle("MEx [GeV]",1);
00263   hMEy                = _dbe->book1D("METTask_MEy",   "METTask_MEy"   ,200,-500,500);
00264   hMEy->setAxisTitle("MEy [GeV]",1);
00265   //hEz                 = _dbe->book1D("METTask_Ez",    "METTask_Ez"    ,500,-500,500);
00266   //hEz->setAxisTitle("MEz [GeV]",1);
00267   hMETSig             = _dbe->book1D("METTask_METSig","METTask_METSig",51,0,51);
00268   hMETSig->setAxisTitle("CaloMETSig",1);
00269   hMET                = _dbe->book1D("METTask_MET",   "METTask_MET"   ,200,0,1000); 
00270   hMET->setAxisTitle("MET [GeV]",1);
00271   hMETPhi             = _dbe->book1D("METTask_METPhi","METTask_METPhi",60,-TMath::Pi(),TMath::Pi()); 
00272   hMETPhi->setAxisTitle("METPhi [rad]",1);
00273   hSumET              = _dbe->book1D("METTask_SumET", "METTask_SumET" ,400,0,2000); 
00274   hSumET->setAxisTitle("SumET [GeV]",1);
00275 
00276   hMET_logx           = _dbe->book1D("METTask_MET_logx",   "METTask_MET_logx"   ,40,-1.,7.);
00277   hMET_logx->setAxisTitle("log(MET) [GeV]",1);
00278   hSumET_logx         = _dbe->book1D("METTask_SumET_logx", "METTask_SumET_logx" ,40,-1.,7.);
00279   hSumET_logx->setAxisTitle("log(SumET) [GeV]",1);
00280 
00281   //hMETIonFeedbck      = _dbe->book1D("METTask_METIonFeedbck", "METTask_METIonFeedbck" ,500,0,1000);
00282   //hMETIonFeedbck->setAxisTitle("MET [GeV]",1);
00283   //hMETHPDNoise        = _dbe->book1D("METTask_METHPDNoise",   "METTask_METHPDNoise"   ,500,0,1000);
00284   //hMETHPDNoise->setAxisTitle("MET [GeV]",1);
00285   //hMETRBXNoise        = _dbe->book1D("METTask_METRBXNoise",   "METTask_METRBXNoise"   ,500,0,1000);
00286   //hMETRBXNoise->setAxisTitle("MET [GeV]",1);
00287 
00288   if (_allhist){
00289     if (bLumiSecPlot){
00290       hMExLS = _dbe->book2D("METTask_MEx_LS","METTask_MEx_LS",200,-200,200,50,0.,500.);
00291       hMExLS->setAxisTitle("MEx [GeV]",1);
00292       hMExLS->setAxisTitle("Lumi Section",2);
00293       hMEyLS = _dbe->book2D("METTask_MEy_LS","METTask_MEy_LS",200,-200,200,50,0.,500.);
00294       hMEyLS->setAxisTitle("MEy [GeV]",1);
00295       hMEyLS->setAxisTitle("Lumi Section",2);
00296     }
00297   }
00298 
00299   if (theMETCollectionLabel.label() == "tcMet" ) {
00300     htrkPt    = _dbe->book1D("METTask_trackPt", "METTask_trackPt", 50, 0, 500);
00301     htrkEta   = _dbe->book1D("METTask_trackEta", "METTask_trackEta", 60, -3.0, 3.0);
00302     htrkNhits = _dbe->book1D("METTask_trackNhits", "METTask_trackNhits", 50, 0, 50);
00303     htrkChi2  = _dbe->book1D("METTask_trackNormalizedChi2", "METTask_trackNormalizedChi2", 20, 0, 20);
00304     htrkD0    = _dbe->book1D("METTask_trackD0", "METTask_trackd0", 50, -1, 1);
00305     helePt    = _dbe->book1D("METTask_electronPt", "METTask_electronPt", 50, 0, 500);
00306     heleEta   = _dbe->book1D("METTask_electronEta", "METTask_electronEta", 60, -3.0, 3.0);
00307     heleHoE   = _dbe->book1D("METTask_electronHoverE", "METTask_electronHoverE", 25, 0, 0.5);
00308     hmuPt     = _dbe->book1D("METTask_muonPt", "METTask_muonPt", 50, 0, 500);
00309     hmuEta    = _dbe->book1D("METTask_muonEta", "METTask_muonEta", 60, -3.0, 3.0);
00310     hmuNhits  = _dbe->book1D("METTask_muonNhits", "METTask_muonNhits", 50, 0, 50);
00311     hmuChi2   = _dbe->book1D("METTask_muonNormalizedChi2", "METTask_muonNormalizedChi2", 20, 0, 20);
00312     hmuD0     = _dbe->book1D("METTask_muonD0", "METTask_muonD0", 50, -1, 1);
00313   }
00314 
00315   hMExCorrection       = _dbe->book1D("METTask_MExCorrection", "METTask_MExCorrection", 100, -500.0,500.0);
00316   hMEyCorrection       = _dbe->book1D("METTask_MEyCorrection", "METTask_MEyCorrection", 100, -500.0,500.0);
00317   hMuonCorrectionFlag  = _dbe->book1D("METTask_CorrectionFlag","METTask_CorrectionFlag", 5, -0.5, 4.5);
00318 
00319 }
00320 
00321 // ***********************************************************
00322 void METAnalyzer::beginRun(const edm::Run& iRun, const edm::EventSetup& iSetup)
00323 {
00324   if ( _HighPtJetEventFlag->on() ) _HighPtJetEventFlag->initRun( iRun, iSetup );
00325   if ( _LowPtJetEventFlag ->on() ) _LowPtJetEventFlag ->initRun( iRun, iSetup );
00326   if ( _MinBiasEventFlag  ->on() ) _MinBiasEventFlag  ->initRun( iRun, iSetup );
00327   if ( _HighMETEventFlag  ->on() ) _HighMETEventFlag  ->initRun( iRun, iSetup );
00328   if ( _LowMETEventFlag   ->on() ) _LowMETEventFlag   ->initRun( iRun, iSetup );
00329   if ( _EleEventFlag      ->on() ) _EleEventFlag      ->initRun( iRun, iSetup );
00330   if ( _MuonEventFlag     ->on() ) _MuonEventFlag     ->initRun( iRun, iSetup );
00331 
00332   if (_HighPtJetEventFlag->on() && _HighPtJetEventFlag->expressionsFromDB(_HighPtJetEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00333     highPtJetExpr_ = _HighPtJetEventFlag->expressionsFromDB(_HighPtJetEventFlag->hltDBKey(), iSetup);
00334   if (_LowPtJetEventFlag->on() && _LowPtJetEventFlag->expressionsFromDB(_LowPtJetEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00335     lowPtJetExpr_  = _LowPtJetEventFlag->expressionsFromDB(_LowPtJetEventFlag->hltDBKey(),   iSetup);
00336   if (_HighMETEventFlag->on() && _HighMETEventFlag->expressionsFromDB(_HighMETEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00337     highMETExpr_   = _HighMETEventFlag->expressionsFromDB(_HighMETEventFlag->hltDBKey(),     iSetup);
00338   if (_LowMETEventFlag->on() && _LowMETEventFlag->expressionsFromDB(_LowMETEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00339     lowMETExpr_    = _LowMETEventFlag->expressionsFromDB(_LowMETEventFlag->hltDBKey(),       iSetup);
00340   if (_MuonEventFlag->on() && _MuonEventFlag->expressionsFromDB(_MuonEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00341     muonExpr_      = _MuonEventFlag->expressionsFromDB(_MuonEventFlag->hltDBKey(),           iSetup);
00342   if (_EleEventFlag->on() && _EleEventFlag->expressionsFromDB(_EleEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00343     elecExpr_      = _EleEventFlag->expressionsFromDB(_EleEventFlag->hltDBKey(),             iSetup);
00344   if (_MinBiasEventFlag->on() && _MinBiasEventFlag->expressionsFromDB(_MinBiasEventFlag->hltDBKey(), iSetup)[0] != "CONFIG_ERROR")
00345     minbiasExpr_   = _MinBiasEventFlag->expressionsFromDB(_MinBiasEventFlag->hltDBKey(),     iSetup);
00346 
00347 }
00348 
00349 // ***********************************************************
00350 void METAnalyzer::endRun(const edm::Run& iRun, const edm::EventSetup& iSetup, DQMStore * dbe)
00351 {
00352   
00353   //
00354   //--- Check the time length of the Run from the lumi section plots
00355 
00356   std::string dirName = _FolderName+_source+"/";
00357   _dbe->setCurrentFolder(dirName);
00358 
00359   TH1F* tlumisec;
00360 
00361   MonitorElement *meLumiSec = _dbe->get("aaa");
00362   meLumiSec = _dbe->get("JetMET/lumisec");
00363 
00364   int totlsec=0;
00365   double totltime=0.;
00366   if ( meLumiSec->getRootObject() ) {
00367     tlumisec = meLumiSec->getTH1F();
00368     for (int i=0; i<500; i++){
00369       if (tlumisec->GetBinContent(i+1)) totlsec++;
00370     }
00371     totltime = double(totlsec*90); // one lumi sec ~ 90 (sec)
00372   }
00373 
00374   if (totltime==0.) totltime=1.; 
00375 
00376   //
00377   //--- Make the integrated plots with rate (Hz)
00378 
00379   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); ic != _FolderNames.end(); ic++)
00380     {
00381 
00382       std::string DirName;
00383       DirName = dirName+*ic;
00384 
00385       makeRatePlot(DirName,totltime);
00386       if ( _HighPtJetEventFlag->on() ) 
00387         makeRatePlot(DirName+"/"+"triggerName_HighJetPt",totltime);
00388       if ( _LowPtJetEventFlag->on() ) 
00389         makeRatePlot(DirName+"/"+"triggerName_LowJetPt",totltime);
00390       if ( _MinBiasEventFlag->on() ) 
00391         makeRatePlot(DirName+"/"+"triggerName_MinBias",totltime);
00392       if ( _HighMETEventFlag->on() ) 
00393         makeRatePlot(DirName+"/"+"triggerName_HighMET",totltime);
00394       if ( _LowMETEventFlag->on() ) 
00395         makeRatePlot(DirName+"/"+"triggerName_LowMET",totltime);
00396       if ( _EleEventFlag->on() ) 
00397         makeRatePlot(DirName+"/"+"triggerName_Ele",totltime);
00398       if ( _MuonEventFlag->on() ) 
00399         makeRatePlot(DirName+"/"+"triggerName_Muon",totltime);
00400     }
00401 }
00402 
00403 // ***********************************************************
00404 void METAnalyzer::makeRatePlot(std::string DirName, double totltime)
00405 {
00406 
00407   _dbe->setCurrentFolder(DirName);
00408   MonitorElement *meMET = _dbe->get(DirName+"/"+"METTask_MET");
00409 
00410   TH1F* tMET;
00411   TH1F* tMETRate;
00412 
00413   if ( meMET )
00414     if ( meMET->getRootObject() ) {
00415       tMET     = meMET->getTH1F();
00416       
00417       // Integral plot & convert number of events to rate (hz)
00418       tMETRate = (TH1F*) tMET->Clone("METTask_METRate");
00419       for (int i = tMETRate->GetNbinsX()-1; i>=0; i--){
00420         tMETRate->SetBinContent(i+1,tMETRate->GetBinContent(i+2)+tMET->GetBinContent(i+1));
00421       }
00422       for (int i = 0; i<tMETRate->GetNbinsX(); i++){
00423         tMETRate->SetBinContent(i+1,tMETRate->GetBinContent(i+1)/double(totltime));
00424       }      
00425 
00426       tMETRate->SetName("METTask_METRate");
00427       tMETRate->SetTitle("METTask_METRate");
00428       hMETRate      = _dbe->book1D("METTask_METRate",tMETRate);
00429     }
00430 }
00431 
00432 // ***********************************************************
00433 void METAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup, 
00434                             const edm::TriggerResults& triggerResults) {
00435 
00436   if (_verbose) std::cout << "METAnalyzer analyze" << std::endl;
00437 
00438   std::string DirName = _FolderName+_source;
00439 
00440   LogTrace(metname)<<"[METAnalyzer] Analyze MET";
00441 
00442   hmetME->Fill(2);
00443 
00444   // ==========================================================  
00445   // Trigger information 
00446   //
00447   _trig_JetMB=0;
00448   _trig_HighPtJet=0;
00449   _trig_LowPtJet=0;
00450   _trig_MinBias=0;
00451   _trig_HighMET=0;
00452   _trig_LowMET=0;
00453   _trig_Ele=0;
00454   _trig_Muon=0;
00455   _trig_PhysDec=0;
00456   if(&triggerResults) {   
00457 
00459 
00460     //
00461     //
00462     // Check how many HLT triggers are in triggerResults 
00463     int ntrigs = triggerResults.size();
00464     if (_verbose) std::cout << "ntrigs=" << ntrigs << std::endl;
00465     
00466     //
00467     //
00468     // If index=ntrigs, this HLT trigger doesn't exist in the HLT table for this data.
00469     const edm::TriggerNames & triggerNames = iEvent.triggerNames(triggerResults);
00470     
00471     //
00472     //
00473     const unsigned int nTrig(triggerNames.size());
00474     for (unsigned int i=0;i<nTrig;++i)
00475       {
00476         if (triggerNames.triggerName(i).find(highPtJetExpr_[0].substr(0,highPtJetExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00477           _trig_HighPtJet=true;
00478         else if (triggerNames.triggerName(i).find(lowPtJetExpr_[0].substr(0,lowPtJetExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00479           _trig_LowPtJet=true;
00480         else if (triggerNames.triggerName(i).find(highMETExpr_[0].substr(0,highMETExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00481           _trig_HighMET=true;
00482         else if (triggerNames.triggerName(i).find(lowMETExpr_[0].substr(0,lowMETExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00483           _trig_LowMET=true;
00484         else if (triggerNames.triggerName(i).find(muonExpr_[0].substr(0,muonExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00485           _trig_Muon=true;
00486         else if (triggerNames.triggerName(i).find(elecExpr_[0].substr(0,elecExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00487           _trig_Ele=true;
00488         else if (triggerNames.triggerName(i).find(minbiasExpr_[0].substr(0,minbiasExpr_[0].rfind("_v")+2))!=std::string::npos && triggerResults.accept(i))
00489           _trig_MinBias=true;
00490       }
00491 
00492     // count number of requested Jet or MB HLT paths which have fired
00493     for (unsigned int i=0; i!=HLTPathsJetMBByName_.size(); i++) {
00494       unsigned int triggerIndex = triggerNames.triggerIndex(HLTPathsJetMBByName_[i]);
00495       if (triggerIndex<triggerResults.size()) {
00496         if (triggerResults.accept(triggerIndex)) {
00497           _trig_JetMB++;
00498         }
00499       }
00500     }
00501     // for empty input vectors (n==0), take all HLT triggers!
00502     if (HLTPathsJetMBByName_.size()==0) _trig_JetMB=triggerResults.size()-1;
00503 
00504     /*
00505       if ( _HighPtJetEventFlag->on() && _HighPtJetEventFlag->accept( iEvent, iSetup) )
00506       _trig_HighPtJet=1;
00507       
00508       if ( _LowPtJetEventFlag->on() && _LowPtJetEventFlag->accept( iEvent, iSetup) )
00509       _trig_LowPtJet=1;
00510       
00511       if ( _MinBiasEventFlag->on() && _MinBiasEventFlag->accept( iEvent, iSetup) )
00512       _trig_MinBias=1;
00513       
00514       if ( _HighMETEventFlag->on() && _HighMETEventFlag->accept( iEvent, iSetup) )
00515       _trig_HighMET=1;
00516       
00517       if ( _LowMETEventFlag->on() && _LowMETEventFlag->accept( iEvent, iSetup) )
00518       _trig_LowMET=1;
00519       
00520       if ( _EleEventFlag->on() && _EleEventFlag->accept( iEvent, iSetup) )
00521       _trig_Ele=1;
00522       
00523       if ( _MuonEventFlag->on() && _MuonEventFlag->accept( iEvent, iSetup) )
00524       _trig_Muon=1;
00525     */
00526     if (triggerNames.triggerIndex(_hlt_PhysDec)   != triggerNames.size() &&
00527         triggerResults.accept(triggerNames.triggerIndex(_hlt_PhysDec)))   _trig_PhysDec=1;
00528   } else {
00529 
00530     edm::LogInfo("MetAnalyzer") << "TriggerResults::HLT not found, "
00531       "automatically select events"; 
00532 
00533     // TriggerResults object not found. Look at all events.    
00534     _trig_JetMB=1;
00535   }
00536 
00537   // ==========================================================
00538   // MET information
00539   
00540   // **** Get the MET container  
00541   edm::Handle<reco::METCollection> metcoll;
00542   iEvent.getByLabel(theMETCollectionLabel, metcoll);
00543   
00544   if(!metcoll.isValid()) {
00545     std::cout<<"Unable to find MET results for MET collection "<<theMETCollectionLabel<<std::endl;
00546     return;
00547   }
00548 
00549   const METCollection *metcol = metcoll.product();
00550   const MET *met;
00551   met = &(metcol->front());
00552     
00553   LogTrace(metname)<<"[METAnalyzer] Call to the MET analyzer";
00554 
00555   // ==========================================================
00556   // TCMET 
00557 
00558   if (theMETCollectionLabel.label() == "tcMet" ) {
00559 
00560     iEvent.getByLabel(inputMuonLabel, muon_h);
00561     iEvent.getByLabel(inputTrackLabel, track_h);
00562     iEvent.getByLabel(inputElectronLabel, electron_h);
00563     iEvent.getByLabel(inputBeamSpotLabel, beamSpot_h);
00564     iEvent.getByLabel("muonTCMETValueMapProducer" , "muCorrData", tcMet_ValueMap_Handle);
00565 
00566     if(!muon_h.isValid())     edm::LogInfo("OutputInfo") << "falied to retrieve muon data require by MET Task";
00567     if(!track_h.isValid())    edm::LogInfo("OutputInfo") << "falied to retrieve track data require by MET Task";
00568     if(!electron_h.isValid()) edm::LogInfo("OutputInfo") << "falied to retrieve electron data require by MET Task";
00569     if(!beamSpot_h.isValid()) edm::LogInfo("OutputInfo") << "falied to retrieve beam spot data require by MET Task";
00570 
00571     bspot = ( beamSpot_h.isValid() ) ? beamSpot_h->position() : math::XYZPoint(0, 0, 0);
00572     
00573   }
00574 
00575   // ==========================================================
00576   // HcalNoiseSummary
00577   //
00578 
00579   edm::Handle<HcalNoiseRBXCollection> HRBXCollection;
00580   iEvent.getByLabel(HcalNoiseRBXCollectionTag,HRBXCollection);
00581   if (!HRBXCollection.isValid()) {
00582     LogDebug("") << "METAnalyzer: Could not find HcalNoiseRBX Collection" << std::endl;
00583     if (_verbose) std::cout << "METAnalyzer: Could not find HcalNoiseRBX Collection" << std::endl;
00584   }
00585   
00586   edm::Handle<HcalNoiseSummary> HNoiseSummary;
00587   iEvent.getByLabel(HcalNoiseSummaryTag,HNoiseSummary);
00588   if (!HNoiseSummary.isValid()) {
00589     LogDebug("") << "METAnalyzer: Could not find Hcal NoiseSummary product" << std::endl;
00590     if (_verbose) std::cout << "METAnalyzer: Could not find Hcal NoiseSummary product" << std::endl;
00591   }
00592 
00593   edm::Handle<reco::CaloJetCollection> caloJets;
00594   iEvent.getByLabel(theJetCollectionLabel, caloJets);
00595   if (!caloJets.isValid()) {
00596     LogDebug("") << "METAnalyzer: Could not find jet product" << std::endl;
00597     if (_verbose) std::cout << "METAnalyzer: Could not find jet product" << std::endl;
00598   }
00599 
00600   // ==========================================================
00601   // MET sanity check
00602 
00603   //   if (_source=="MET") validateMET(*met, tcCandidates);
00604   
00605   // ==========================================================
00606   // JetID 
00607 
00608   if (_verbose) std::cout << "JetID starts" << std::endl;
00609   
00610   //
00611   // --- Minimal cuts
00612   //
00613   bool bJetIDMinimal=true;
00614   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00615        cal!=caloJets->end(); ++cal){
00616     jetID->calculate(iEvent, *cal);
00617     if (cal->pt()>10.){
00618       if (fabs(cal->eta())<=2.6 && 
00619           cal->emEnergyFraction()<=0.01) bJetIDMinimal=false;
00620     }
00621   }
00622 
00623   //
00624   // --- Loose cuts, not  specific for now!
00625   //
00626   bool bJetIDLoose=true;
00627   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00628        cal!=caloJets->end(); ++cal){ 
00629     jetID->calculate(iEvent, *cal);
00630     if (_verbose) std::cout << jetID->n90Hits() << " " 
00631                             << jetID->restrictedEMF() << " "
00632                             << cal->pt() << std::endl;
00633     if (cal->pt()>10.){
00634       //
00635       // for all regions
00636       if (jetID->n90Hits()<2)  bJetIDLoose=false; 
00637       if (jetID->fHPD()>=0.98) bJetIDLoose=false; 
00638       //if (jetID->restrictedEMF()<0.01) bJetIDLoose=false; 
00639       //
00640       // for non-forward
00641       if (fabs(cal->eta())<2.55){
00642         if (cal->emEnergyFraction()<=0.01) bJetIDLoose=false; 
00643       }
00644       // for forward
00645       else {
00646         if (cal->emEnergyFraction()<=-0.9) bJetIDLoose=false; 
00647         if (cal->pt()>80.){
00648         if (cal->emEnergyFraction()>= 1.0) bJetIDLoose=false; 
00649         }
00650       } // forward vs non-forward
00651     }   // pt>10 GeV/c
00652   }     // calor-jets loop
00653  
00654   //
00655   // --- Tight cuts
00656   //
00657   bool bJetIDTight=true;
00658   bJetIDTight=bJetIDLoose;
00659   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
00660        cal!=caloJets->end(); ++cal){
00661     jetID->calculate(iEvent, *cal);
00662     if (cal->pt()>25.){
00663       //
00664       // for all regions
00665       if (jetID->fHPD()>=0.95) bJetIDTight=false; 
00666       //
00667       // for 1.0<|eta|<1.75
00668       if (fabs(cal->eta())>=1.00 && fabs(cal->eta())<1.75){
00669         if (cal->pt()>80. && cal->emEnergyFraction()>=1.) bJetIDTight=false; 
00670       }
00671       //
00672       // for 1.75<|eta|<2.55
00673       else if (fabs(cal->eta())>=1.75 && fabs(cal->eta())<2.55){
00674         if (cal->pt()>80. && cal->emEnergyFraction()>=1.) bJetIDTight=false; 
00675       }
00676       //
00677       // for 2.55<|eta|<3.25
00678       else if (fabs(cal->eta())>=2.55 && fabs(cal->eta())<3.25){
00679         if (cal->pt()< 50.                   && cal->emEnergyFraction()<=-0.3) bJetIDTight=false; 
00680         if (cal->pt()>=50. && cal->pt()< 80. && cal->emEnergyFraction()<=-0.2) bJetIDTight=false; 
00681         if (cal->pt()>=80. && cal->pt()<340. && cal->emEnergyFraction()<=-0.1) bJetIDTight=false; 
00682         if (cal->pt()>=340.                  && cal->emEnergyFraction()<=-0.1 
00683                                              && cal->emEnergyFraction()>=0.95) bJetIDTight=false; 
00684       }
00685       //
00686       // for 3.25<|eta|
00687       else if (fabs(cal->eta())>=3.25){
00688         if (cal->pt()< 50.                   && cal->emEnergyFraction()<=-0.3
00689                                              && cal->emEnergyFraction()>=0.90) bJetIDTight=false; 
00690         if (cal->pt()>=50. && cal->pt()<130. && cal->emEnergyFraction()<=-0.2
00691                                              && cal->emEnergyFraction()>=0.80) bJetIDTight=false; 
00692         if (cal->pt()>=130.                  && cal->emEnergyFraction()<=-0.1 
00693                                              && cal->emEnergyFraction()>=0.70) bJetIDTight=false; 
00694       }
00695     }   // pt>10 GeV/c
00696   }     // calor-jets loop
00697   
00698   if (_verbose) std::cout << "JetID ends" << std::endl;
00699 
00700 
00701   // ==========================================================
00702   // HCAL Noise filter
00703   
00704   bool bHcalNoiseFilter      = HNoiseSummary->passLooseNoiseFilter();
00705   bool bHcalNoiseFilterTight = HNoiseSummary->passTightNoiseFilter();
00706 
00707   // ==========================================================
00708   // Get BeamHaloSummary
00709   edm::Handle<BeamHaloSummary> TheBeamHaloSummary ;
00710   iEvent.getByLabel(BeamHaloSummaryTag, TheBeamHaloSummary) ;
00711 
00712   if (!TheBeamHaloSummary.isValid()) {
00713     std::cout << "BeamHaloSummary doesn't exist" << std::endl;
00714   }
00715 
00716   bool bBeamHaloIDTightPass = true;
00717   bool bBeamHaloIDLoosePass = true;
00718 
00719   if(!TheBeamHaloSummary.isValid()) {
00720 
00721   const BeamHaloSummary TheSummary = (*TheBeamHaloSummary.product() );
00722 
00723   if( !TheSummary.EcalLooseHaloId()  && !TheSummary.HcalLooseHaloId() && 
00724       !TheSummary.CSCLooseHaloId()   && !TheSummary.GlobalLooseHaloId() )
00725     bBeamHaloIDLoosePass = false;
00726 
00727   if( !TheSummary.EcalTightHaloId()  && !TheSummary.HcalTightHaloId() && 
00728       !TheSummary.CSCTightHaloId()   && !TheSummary.GlobalTightHaloId() )
00729     bBeamHaloIDTightPass = false;
00730 
00731   }
00732 
00733   // ==========================================================
00734   //Vertex information
00735 
00736   bool bPrimaryVertex = true;
00737   if(_doPVCheck){
00738     bPrimaryVertex = false;
00739     Handle<VertexCollection> vertexHandle;
00740 
00741     iEvent.getByLabel(vertexTag, vertexHandle);
00742 
00743     if (!vertexHandle.isValid()) {
00744       LogDebug("") << "CaloMETAnalyzer: Could not find vertex collection" << std::endl;
00745       if (_verbose) std::cout << "CaloMETAnalyzer: Could not find vertex collection" << std::endl;
00746     }
00747     
00748     if ( vertexHandle.isValid() ){
00749       VertexCollection vertexCollection = *(vertexHandle.product());
00750       int vertex_number     = vertexCollection.size();
00751       VertexCollection::const_iterator v = vertexCollection.begin();
00752       for ( ; v != vertexCollection.end(); ++v) {
00753         double vertex_chi2    = v->normalizedChi2();
00754         double vertex_ndof    = v->ndof();
00755         bool   fakeVtx        = v->isFake();
00756         double vertex_Z       = v->z();
00757         
00758         if (  !fakeVtx
00759               && vertex_number>=_nvtx_min
00760               && vertex_ndof   >_vtxndof_min
00761               && vertex_chi2   <_vtxchi2_max
00762               && fabs(vertex_Z)<_vtxz_max ) 
00763           bPrimaryVertex = true;
00764       }
00765     }
00766   }
00767   // ==========================================================
00768 
00769   edm::Handle< L1GlobalTriggerReadoutRecord > gtReadoutRecord;
00770   iEvent.getByLabel( gtTag, gtReadoutRecord);
00771 
00772   if (!gtReadoutRecord.isValid()) {
00773     LogDebug("") << "CaloMETAnalyzer: Could not find GT readout record" << std::endl;
00774     if (_verbose) std::cout << "CaloMETAnalyzer: Could not find GT readout record product" << std::endl;
00775   }
00776   
00777   bool bTechTriggers    = true;
00778   bool bTechTriggersAND = true;
00779   bool bTechTriggersOR  = false;
00780   bool bTechTriggersNOT = false;
00781 
00782   if (gtReadoutRecord.isValid()) {
00783     const TechnicalTriggerWord&  technicalTriggerWordBeforeMask = gtReadoutRecord->technicalTriggerWord();
00784     
00785     if (_techTrigsAND.size() == 0)
00786       bTechTriggersAND = true;
00787     else
00788       for (unsigned ttr = 0; ttr != _techTrigsAND.size(); ttr++) {
00789         bTechTriggersAND = bTechTriggersAND && technicalTriggerWordBeforeMask.at(_techTrigsAND.at(ttr));
00790       }
00791     
00792     if (_techTrigsAND.size() == 0)
00793       bTechTriggersOR = true;
00794     else
00795       for (unsigned ttr = 0; ttr != _techTrigsOR.size(); ttr++) {
00796         bTechTriggersOR = bTechTriggersOR || technicalTriggerWordBeforeMask.at(_techTrigsOR.at(ttr));
00797       }
00798     if (_techTrigsNOT.size() == 0)
00799       bTechTriggersNOT = false;
00800     else
00801       for (unsigned ttr = 0; ttr != _techTrigsNOT.size(); ttr++) {
00802         bTechTriggersNOT = bTechTriggersNOT || technicalTriggerWordBeforeMask.at(_techTrigsNOT.at(ttr));
00803       }
00804   }
00805   else
00806     {
00807       bTechTriggersAND = true;
00808       bTechTriggersOR  = true;
00809       bTechTriggersNOT = false;
00810     }
00811   
00812   if (_techTrigsAND.size()==0)
00813     bTechTriggersAND = true;
00814   if (_techTrigsOR.size()==0)
00815     bTechTriggersOR = true;
00816   if (_techTrigsNOT.size()==0)
00817     bTechTriggersNOT = false;
00818   
00819   bTechTriggers = bTechTriggersAND && bTechTriggersOR && !bTechTriggersNOT;
00820   
00821   // ==========================================================
00822   // Reconstructed MET Information - fill MonitorElements
00823   
00824   bool bHcalNoise   = bHcalNoiseFilter;
00825   bool bBeamHaloID  = bBeamHaloIDLoosePass;
00826   bool bJetID       = true;
00827 
00828   bool bPhysicsDeclared = true;
00829   if(_doHLTPhysicsOn) bPhysicsDeclared =_trig_PhysDec;
00830 
00831 
00832   if      (_tightHcalFiltering)     bHcalNoise  = bHcalNoiseFilterTight;
00833   if      (_tightBHFiltering)       bBeamHaloID = bBeamHaloIDTightPass;
00834 
00835   if      (_tightJetIDFiltering==1)  bJetID      = bJetIDMinimal;
00836   else if (_tightJetIDFiltering==2)  bJetID      = bJetIDLoose;
00837   else if (_tightJetIDFiltering==3)  bJetID      = bJetIDTight;
00838   else if (_tightJetIDFiltering==-1) bJetID      = true;
00839 
00840   bool bBasicCleanup = bTechTriggers && bPrimaryVertex && bPhysicsDeclared;
00841   bool bExtraCleanup = bBasicCleanup && bHcalNoise && bJetID && bBeamHaloID;
00842 
00843   //std::string DirName = _FolderName+_source;
00844   
00845   for (std::vector<std::string>::const_iterator ic = _FolderNames.begin(); 
00846        ic != _FolderNames.end(); ic++){
00847     if (*ic=="All")                                             fillMESet(iEvent, DirName+"/"+*ic, *met);
00848     if (DCSFilter->filter(iEvent, iSetup)) {
00849     if (_cleanupSelection){
00850     if (*ic=="BasicCleanup" && bBasicCleanup)                   fillMESet(iEvent, DirName+"/"+*ic, *met);
00851     if (*ic=="ExtraCleanup" && bExtraCleanup)                   fillMESet(iEvent, DirName+"/"+*ic, *met);
00852     }
00853     if (_allSelection) {
00854       if (*ic=="HcalNoiseFilter"      && bHcalNoiseFilter )       fillMESet(iEvent, DirName+"/"+*ic, *met);
00855       if (*ic=="HcalNoiseFilterTight" && bHcalNoiseFilterTight )  fillMESet(iEvent, DirName+"/"+*ic, *met);
00856       if (*ic=="JetIDMinimal"         && bJetIDMinimal)           fillMESet(iEvent, DirName+"/"+*ic, *met);
00857       if (*ic=="JetIDLoose"           && bJetIDLoose)             fillMESet(iEvent, DirName+"/"+*ic, *met);
00858       if (*ic=="JetIDTight"           && bJetIDTight)             fillMESet(iEvent, DirName+"/"+*ic, *met);
00859       if (*ic=="BeamHaloIDTightPass"  && bBeamHaloIDTightPass)    fillMESet(iEvent, DirName+"/"+*ic, *met);
00860       if (*ic=="BeamHaloIDLoosePass"  && bBeamHaloIDLoosePass)    fillMESet(iEvent, DirName+"/"+*ic, *met);
00861       if (*ic=="Triggers"             && bTechTriggers)           fillMESet(iEvent, DirName+"/"+*ic, *met);
00862       if (*ic=="PV"                   && bPrimaryVertex)          fillMESet(iEvent, DirName+"/"+*ic, *met);
00863     }
00864     } // DCS
00865   }
00866 }
00867 
00868 
00869 // ***********************************************************
00870 void METAnalyzer::fillMESet(const edm::Event& iEvent, std::string DirName, 
00871                               const reco::MET& met)
00872 {
00873 
00874   _dbe->setCurrentFolder(DirName);
00875 
00876   bool bLumiSecPlot=false;
00877   if (DirName.find("All")) bLumiSecPlot=true;
00878 
00879   if (_trig_JetMB)
00880     fillMonitorElement(iEvent,DirName,"",met, bLumiSecPlot);
00881   if (_trig_HighPtJet)
00882     fillMonitorElement(iEvent,DirName,"HighPtJet",met,false);
00883   if (_trig_LowPtJet)
00884     fillMonitorElement(iEvent,DirName,"LowPtJet",met,false);
00885   if (_trig_MinBias)
00886     fillMonitorElement(iEvent,DirName,"MinBias",met,false);
00887   if (_trig_HighMET)
00888     fillMonitorElement(iEvent,DirName,"HighMET",met,false);
00889   if (_trig_LowMET)
00890     fillMonitorElement(iEvent,DirName,"LowMET",met,false);
00891   if (_trig_Ele)
00892     fillMonitorElement(iEvent,DirName,"Ele",met,false);
00893   if (_trig_Muon)
00894     fillMonitorElement(iEvent,DirName,"Muon",met,false);
00895 }
00896 
00897 // ***********************************************************
00898 void METAnalyzer::fillMonitorElement(const edm::Event& iEvent, std::string DirName, 
00899                                          std::string TriggerTypeName, 
00900                                          const reco::MET& met, bool bLumiSecPlot)
00901 {
00902 
00903   if (TriggerTypeName=="HighPtJet") {
00904     if (!selectHighPtJetEvent(iEvent)) return;
00905   }
00906   else if (TriggerTypeName=="LowPtJet") {
00907     if (!selectLowPtJetEvent(iEvent)) return;
00908   }
00909   else if (TriggerTypeName=="HighMET") {
00910     if (met.pt()<_highMETThreshold) return;
00911   }
00912   else if (TriggerTypeName=="LowMET") {
00913     if (met.pt()<_lowMETThreshold) return;
00914   }
00915   else if (TriggerTypeName=="Ele") {
00916     if (!selectWElectronEvent(iEvent)) return;
00917   }
00918   else if (TriggerTypeName=="Muon") {
00919     if (!selectWMuonEvent(iEvent)) return;
00920   }
00921   
00922 // Reconstructed MET Information
00923   double SumET  = met.sumEt();
00924   double METSig = met.mEtSig();
00925   //double Ez     = met.e_longitudinal();
00926   double MET    = met.pt();
00927   double MEx    = met.px();
00928   double MEy    = met.py();
00929   double METPhi = met.phi();
00930 
00931   //
00932   int myLuminosityBlock;
00933   //  myLuminosityBlock = (evtCounter++)/1000;
00934   myLuminosityBlock = iEvent.luminosityBlock();
00935   //
00936 
00937   if (TriggerTypeName!="") DirName = DirName +"/"+TriggerTypeName;
00938 
00939   if (_verbose) std::cout << "_etThreshold = " << _etThreshold << std::endl;
00940   if (SumET>_etThreshold){
00941     
00942     hMEx    = _dbe->get(DirName+"/"+"METTask_MEx");     if (hMEx           && hMEx->getRootObject())     hMEx          ->Fill(MEx);
00943     hMEy    = _dbe->get(DirName+"/"+"METTask_MEy");     if (hMEy           && hMEy->getRootObject())     hMEy          ->Fill(MEy);
00944     hMET    = _dbe->get(DirName+"/"+"METTask_MET");     if (hMET           && hMET->getRootObject())     hMET          ->Fill(MET);
00945     hMETPhi = _dbe->get(DirName+"/"+"METTask_METPhi");  if (hMETPhi        && hMETPhi->getRootObject())  hMETPhi       ->Fill(METPhi);
00946     hSumET  = _dbe->get(DirName+"/"+"METTask_SumET");   if (hSumET         && hSumET->getRootObject())   hSumET        ->Fill(SumET);
00947     hMETSig = _dbe->get(DirName+"/"+"METTask_METSig");  if (hMETSig        && hMETSig->getRootObject())  hMETSig       ->Fill(METSig);
00948     //hEz     = _dbe->get(DirName+"/"+"METTask_Ez");      if (hEz            && hEz->getRootObject())      hEz           ->Fill(Ez);
00949 
00950     hMET_logx   = _dbe->get(DirName+"/"+"METTask_MET_logx");    if (hMET_logx      && hMET_logx->getRootObject())    hMET_logx->Fill(log10(MET));
00951     hSumET_logx = _dbe->get(DirName+"/"+"METTask_SumET_logx");  if (hSumET_logx    && hSumET_logx->getRootObject())  hSumET_logx->Fill(log10(SumET));
00952 
00953     //hMETIonFeedbck = _dbe->get(DirName+"/"+"METTask_METIonFeedbck");  if (hMETIonFeedbck && hMETIonFeedbck->getRootObject())  hMETIonFeedbck->Fill(MET);
00954     //hMETHPDNoise   = _dbe->get(DirName+"/"+"METTask_METHPDNoise");    if (hMETHPDNoise   && hMETHPDNoise->getRootObject())    hMETHPDNoise->Fill(MET);
00955        
00956     if (_allhist){
00957       if (bLumiSecPlot){
00958         hMExLS = _dbe->get(DirName+"/"+"METTask_MExLS"); if (hMExLS  &&  hMExLS->getRootObject())   hMExLS->Fill(MEx,myLuminosityBlock);
00959         hMEyLS = _dbe->get(DirName+"/"+"METTask_MEyLS"); if (hMEyLS  &&  hMEyLS->getRootObject())   hMEyLS->Fill(MEy,myLuminosityBlock);
00960       }
00961     } // _allhist
00962 
00964     if (theMETCollectionLabel.label() == "tcMet" ) {
00965     
00966       if(track_h.isValid()) {
00967         for( edm::View<reco::Track>::const_iterator trkit = track_h->begin(); trkit != track_h->end(); trkit++ ) {
00968           htrkPt    = _dbe->get(DirName+"/"+"METTask_trackPt");     if (htrkPt    && htrkPt->getRootObject())     htrkPt->Fill( trkit->pt() );
00969           htrkEta   = _dbe->get(DirName+"/"+"METTask_trackEta");    if (htrkEta   && htrkEta->getRootObject())    htrkEta->Fill( trkit->eta() );
00970           htrkNhits = _dbe->get(DirName+"/"+"METTask_trackNhits");  if (htrkNhits && htrkNhits->getRootObject())  htrkNhits->Fill( trkit->numberOfValidHits() );
00971           htrkChi2  = _dbe->get(DirName+"/"+"METTask_trackNormalizedChi2");   if (htrkChi2  && htrkChi2->getRootObject())   htrkChi2->Fill( trkit->chi2() / trkit->ndof() );
00972           double d0 = -1 * trkit->dxy( bspot );
00973           htrkD0    = _dbe->get(DirName+"/"+"METTask_trackD0");     if (htrkD0 && htrkD0->getRootObject())        htrkD0->Fill( d0 );
00974         }
00975       }
00976       
00977       if(electron_h.isValid()) {
00978         for( edm::View<reco::GsfElectron>::const_iterator eleit = electron_h->begin(); eleit != electron_h->end(); eleit++ ) {
00979           helePt  = _dbe->get(DirName+"/"+"METTask_electronPt");   if (helePt  && helePt->getRootObject())   helePt->Fill( eleit->p4().pt() );  
00980           heleEta = _dbe->get(DirName+"/"+"METTask_electronEta");  if (heleEta && heleEta->getRootObject())  heleEta->Fill( eleit->p4().eta() );
00981           heleHoE = _dbe->get(DirName+"/"+"METTask_electronHoverE");  if (heleHoE && heleHoE->getRootObject())  heleHoE->Fill( eleit->hadronicOverEm() );
00982         }
00983       }
00984       
00985       if(muon_h.isValid()) {
00986         for( reco::MuonCollection::const_iterator muonit = muon_h->begin(); muonit != muon_h->end(); muonit++ ) {      
00987           const reco::TrackRef siTrack = muonit->innerTrack();      
00988           hmuPt    = _dbe->get(DirName+"/"+"METTask_muonPt");     if (hmuPt    && hmuPt->getRootObject())  hmuPt   ->Fill( muonit->p4().pt() );
00989           hmuEta   = _dbe->get(DirName+"/"+"METTask_muonEta");    if (hmuEta   && hmuEta->getRootObject())  hmuEta  ->Fill( muonit->p4().eta() );
00990           hmuNhits = _dbe->get(DirName+"/"+"METTask_muonNhits");  if (hmuNhits && hmuNhits->getRootObject())  hmuNhits->Fill( siTrack.isNonnull() ? siTrack->numberOfValidHits() : -999 );
00991           hmuChi2  = _dbe->get(DirName+"/"+"METTask_muonNormalizedChi2");   if (hmuChi2  && hmuChi2->getRootObject())  hmuChi2 ->Fill( siTrack.isNonnull() ? siTrack->chi2()/siTrack->ndof() : -999 );
00992           double d0 = siTrack.isNonnull() ? -1 * siTrack->dxy( bspot) : -999;
00993           hmuD0    = _dbe->get(DirName+"/"+"METTask_muonD0");     if (hmuD0    && hmuD0->getRootObject())  hmuD0->Fill( d0 );
00994         }
00995         
00996         const unsigned int nMuons = muon_h->size();      
00997         for( unsigned int mus = 0; mus < nMuons; mus++ ) {
00998           reco::MuonRef muref( muon_h, mus);
00999           reco::MuonMETCorrectionData muCorrData = (*tcMet_ValueMap_Handle)[muref];
01000           hMExCorrection      = _dbe->get(DirName+"/"+"METTask_MExCorrection");       if (hMExCorrection      && hMExCorrection->getRootObject())       hMExCorrection-> Fill(muCorrData.corrY());
01001           hMEyCorrection      = _dbe->get(DirName+"/"+"METTask_MEyCorrection");       if (hMEyCorrection      && hMEyCorrection->getRootObject())       hMEyCorrection-> Fill(muCorrData.corrX());
01002           hMuonCorrectionFlag = _dbe->get(DirName+"/"+"METTask_CorrectionFlag");  if (hMuonCorrectionFlag && hMuonCorrectionFlag->getRootObject())  hMuonCorrectionFlag-> Fill(muCorrData.type());
01003         }
01004       }
01005     }
01006 
01007   } // et threshold cut
01008 
01009 }
01010 
01011 // ***********************************************************
01012 bool METAnalyzer::selectHighPtJetEvent(const edm::Event& iEvent){
01013 
01014   bool return_value=false;
01015 
01016   edm::Handle<reco::CaloJetCollection> caloJets;
01017   iEvent.getByLabel(theJetCollectionLabel, caloJets);
01018   if (!caloJets.isValid()) {
01019     LogDebug("") << "METAnalyzer: Could not find jet product" << std::endl;
01020     if (_verbose) std::cout << "METAnalyzer: Could not find jet product" << std::endl;
01021   }
01022 
01023   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
01024        cal!=caloJets->end(); ++cal){
01025     if (cal->pt()>_highPtJetThreshold){
01026       return_value=true;
01027     }
01028   }
01029   
01030   return return_value;
01031 }
01032 
01033 // // ***********************************************************
01034 bool METAnalyzer::selectLowPtJetEvent(const edm::Event& iEvent){
01035 
01036   bool return_value=false;
01037 
01038   edm::Handle<reco::CaloJetCollection> caloJets;
01039   iEvent.getByLabel(theJetCollectionLabel, caloJets);
01040   if (!caloJets.isValid()) {
01041     LogDebug("") << "METAnalyzer: Could not find jet product" << std::endl;
01042     if (_verbose) std::cout << "METAnalyzer: Could not find jet product" << std::endl;
01043   }
01044 
01045   for (reco::CaloJetCollection::const_iterator cal = caloJets->begin(); 
01046        cal!=caloJets->end(); ++cal){
01047     if (cal->pt()>_lowPtJetThreshold){
01048       return_value=true;
01049     }
01050   }
01051 
01052   return return_value;
01053 
01054 }
01055 
01056 
01057 // ***********************************************************
01058 bool METAnalyzer::selectWElectronEvent(const edm::Event& iEvent){
01059 
01060   bool return_value=true;
01061 
01062   /*
01063     W-electron event selection comes here
01064    */
01065 
01066   return return_value;
01067 
01068 }
01069 
01070 // ***********************************************************
01071 bool METAnalyzer::selectWMuonEvent(const edm::Event& iEvent){
01072 
01073   bool return_value=true;
01074 
01075   /*
01076     W-muon event selection comes here
01077    */
01078 
01079   return return_value;
01080 
01081 }
01082
/data/refman/pasoursint/CMSSW_4_4_5_patch3/src/DQMOffline/JetMET/src/METAnalyzer.cc
