Reference Manual

Go to the documentation of this file.
00001 #include "DQM/Physics/src/HiggsDQM.h"
00002 
00003 #include <memory>
00004 
00005 // Framework
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 #include "FWCore/Framework/interface/EventSetup.h"
00008 #include "FWCore/Framework/interface/LuminosityBlock.h"
00009 #include "FWCore/ServiceRegistry/interface/Service.h"
00010 #include "FWCore/ParameterSet/interface/FileInPath.h"
00011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00012 #include "FWCore/Framework/interface/Run.h"
00013 #include "DataFormats/Provenance/interface/EventID.h"
00014 
00015 // Candidate handling
00016 #include "DataFormats/Candidate/interface/Candidate.h"
00017 #include "DataFormats/Candidate/interface/CandidateFwd.h"
00018 #include "DataFormats/Candidate/interface/OverlapChecker.h"
00019 #include "DataFormats/Candidate/interface/CompositeCandidate.h"
00020 #include "DataFormats/Candidate/interface/CompositeCandidateFwd.h"
00021 #include "DataFormats/Candidate/interface/CandMatchMap.h"
00022 #include <DataFormats/EgammaCandidates/interface/GsfElectron.h>
00023 #include <DataFormats/EgammaCandidates/interface/GsfElectronFwd.h>
00024 #include "DataFormats/MuonReco/interface/Muon.h"
00025 #include "DataFormats/MuonReco/interface/MuonFwd.h"
00026 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
00027 
00028 #include "DQMServices/Core/interface/DQMStore.h"
00029 #include "DQMServices/Core/interface/MonitorElement.h"
00030 
00031 //#include "HiggsAnalysis/HiggsToZZ4Leptons/plugins/HZZ4LeptonsElectronAssociationMap.h"
00032 //#include "HiggsAnalysis/HiggsToZZ4Leptons/plugins/HZZ4LeptonsMuonAssociationMap.h"
00033 
00034 // vertexing
00035 // Transient tracks
00036 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00037 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
00038 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
00039 #include "RecoTracker/Record/interface/TrackerRecoGeometryRecord.h"
00040 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00041 #include "Geometry/Records/interface/MuonGeometryRecord.h"
00042 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00043 
00044 // Vertex utilities
00045 #include "DataFormats/VertexReco/interface/Vertex.h"
00046 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00047 
00048 // Geometry
00049 #include "DataFormats/Math/interface/Vector3D.h"
00050 #include "DataFormats/Math/interface/LorentzVector.h"
00051 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
00052 #include "DataFormats/Common/interface/AssociationVector.h"
00053 #include "DataFormats/Math/interface/Point3D.h"
00054 #include "DataFormats/Math/interface/deltaR.h"
00055 #include "DataFormats/Math/interface/deltaPhi.h"
00056 
00057 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
00058 #include "DataFormats/EgammaCandidates/interface/Electron.h"
00059 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
00060 #include "DataFormats/DetId/interface/DetId.h"
00061 #include "DataFormats/Common/interface/RefToBase.h"
00062 
00063 
00064 // RECO
00065 // Isolation
00066 //#include "HiggsAnalysis/HiggsToZZ4Leptons/plugins/CandidateHadIsolation.h"
00067 //#include "HiggsAnalysis/HiggsToZZ4Leptons/plugins/CandidateTkIsolation.h"
00068 #include "DataFormats/MuonReco/interface/Muon.h"
00069 #include "DataFormats/MuonReco/interface/MuonFwd.h"
00070 #include "DataFormats/MuonReco/interface/MuonIsolation.h" 
00071 #include "DataFormats/TrackReco/interface/Track.h"
00072 #include <DataFormats/EgammaCandidates/interface/GsfElectron.h>
00073 #include "DataFormats/EgammaCandidates/interface/Photon.h"
00074 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
00075 #include "DataFormats/VertexReco/interface/Vertex.h"
00076 #include "DataFormats/JetReco/interface/CaloJet.h"
00077 
00078 //MET
00079 #include "DataFormats/METReco/interface/CaloMETCollection.h"
00080 #include "DataFormats/METReco/interface/PFMETCollection.h"
00081 #include "DataFormats/METReco/interface/GenMETCollection.h"
00082 #include "DataFormats/METReco/interface/CaloMET.h"
00083 #include "DataFormats/METReco/interface/PFMET.h"
00084 #include "DataFormats/METReco/interface/MET.h"
00085 #include "DataFormats/METReco/interface/METCollection.h"
00086 
00087 #include "DataFormats/Math/interface/LorentzVector.h"
00088 #include "TLorentzVector.h"
00089 
00090 #include <iostream>
00091 #include <iomanip>
00092 #include <stdio.h>
00093 #include <string>
00094 #include <sstream>
00095 #include <math.h>
00096 
00097 using namespace edm;
00098 using namespace std;
00099 using namespace reco;
00100 
00101 struct SortCandByDecreasingPt {
00102   bool operator()( const Candidate &c1, const Candidate &c2) const {
00103     return c1.pt() > c2.pt();
00104   }
00105 };
00106 
00107 
00108 double HiggsDQM::Distance( const reco::Candidate & c1, const reco::Candidate & c2 ) {
00109         return  deltaR(c1,c2);
00110 }
00111 
00112 double HiggsDQM::DistancePhi( const reco::Candidate & c1, const reco::Candidate & c2 ) {
00113         return  deltaPhi(c1.p4().phi(),c2.p4().phi());
00114 }
00115 
00116 // This always returns only a positive deltaPhi
00117 double HiggsDQM::calcDeltaPhi(double phi1, double phi2) {
00118   double deltaPhi = phi1 - phi2;
00119   if (deltaPhi < 0) deltaPhi = -deltaPhi;
00120   if (deltaPhi > 3.1415926) {
00121     deltaPhi = 2 * 3.1415926 - deltaPhi;
00122   }
00123   return deltaPhi;
00124 }
00125 
00126 //
00127 // -- Constructor
00128 //
00129 HiggsDQM::HiggsDQM(const edm::ParameterSet& ps){
00130  //cout<<"Entering  HiggsDQM::HiggsDQM: "<<endl;
00131  
00132   edm::LogInfo("HZZ4LeptonsDQM") <<  " Creating HZZ4LeptonsDQM " << "\n" ;
00133 
00134   bei_ = Service<DQMStore>().operator->();
00135   bei_->setCurrentFolder("Physics/Higgs");
00136   bookHistos(bei_);
00137 
00138   typedef std::vector<edm::InputTag> vtag;
00139   // Get parameters from configuration file
00140   theElecTriggerPathToPass    = ps.getParameter<string>("elecTriggerPathToPass");
00141   theMuonTriggerPathToPass    = ps.getParameter<string>("muonTriggerPathToPass");
00142   theTriggerResultsCollection = ps.getParameter<InputTag>("triggerResultsCollection");
00143   theMuonCollectionLabel      = ps.getParameter<InputTag>("muonCollection");
00144   theElectronCollectionLabel  = ps.getParameter<InputTag>("electronCollection");
00145   theCaloJetCollectionLabel   = ps.getParameter<InputTag>("caloJetCollection");
00146   theCaloMETCollectionLabel   = ps.getParameter<InputTag>("caloMETCollection");
00147   thePfMETCollectionLabel     = ps.getParameter<InputTag>("pfMETCollection");
00148   // just to initialize
00149   isValidHltConfig_ = false;
00150   // cuts:
00151   ptThrMu1_ = ps.getUntrackedParameter<double>("PtThrMu1");
00152   ptThrMu2_ = ps.getUntrackedParameter<double>("PtThrMu2");
00153  
00154 
00155   
00156  //cout<<"...leaving  HiggsDQM::HiggsDQM. "<<endl;
00157 }
00158 //
00159 // -- Destructor
00160 //
00161 HiggsDQM::~HiggsDQM(){
00162   //cout<<"Entering HiggsDQM::~HiggsDQM: "<<endl;
00163   
00164   edm::LogInfo("HiggsDQM") <<  " Deleting HiggsDQM " << "\n" ;
00165 
00166   //cout<<"...leaving HiggsDQM::~HiggsDQM. "<<endl;
00167 }
00168 //
00169 // -- Begin Job
00170 //
00171 void HiggsDQM::beginJob(){
00172   //cout<<"Entering HiggsDQM::beginJob: "<<endl;
00173 
00174   nLumiSecs_ = 0;
00175   nEvents_   = 0;
00176   pi = 3.14159265;
00177   
00178 
00179   //cout<<"...leaving HiggsDQM::beginJob. "<<endl;
00180 }
00181 //
00182 // -- Begin Run
00183 //
00184 void HiggsDQM::beginRun(Run const& run, edm::EventSetup const& eSetup) {
00185   edm::LogInfo ("HiggsDQM") <<"[HiggsDQM]: Begining of Run";
00186   // passed as parameter to HLTConfigProvider::init(), not yet used
00187   bool isConfigChanged = false;
00188   
00189   // isValidHltConfig_ used to short-circuit analyze() in case of problems
00190   //  const std::string hltProcessName( "HLT" );
00191   const std::string hltProcessName = theTriggerResultsCollection.process();
00192   isValidHltConfig_ = hltConfigProvider_.init( run, eSetup, hltProcessName, isConfigChanged );
00193 
00194 }
00195 //
00196 // -- Begin  Luminosity Block
00197 //
00198 void HiggsDQM::beginLuminosityBlock(edm::LuminosityBlock const& lumiSeg, 
00199                                             edm::EventSetup const& context) {
00200   //cout<<"Entering HiggsDQM::beginLuminosityBlock: "<<endl;
00201   
00202   edm::LogInfo ("HiggsDQM") <<"[HiggsDQM]: Begin of LS transition";
00203 
00204   //cout<<"...leaving HiggsDQM::beginLuminosityBlock. "<<endl;
00205 }
00206 //
00207 //  -- Book histograms
00208 //
00209 void HiggsDQM::bookHistos(DQMStore* bei){
00210   bei->cd();
00211   bei->setCurrentFolder("Physics/Higgs");
00212   h_vertex_number = bei->book1D("h_vertex_number", "Number of event vertices in collection", 10,-0.5,   9.5 );
00213   h_vertex_chi2  = bei->book1D("h_vertex_chi2" , "Event Vertex #chi^{2}/n.d.o.f."          , 100, 0.0,   2.0 );
00214   h_vertex_numTrks = bei->book1D("h_vertex_numTrks", "Event Vertex, number of tracks"     , 100, -0.5,  99.5 );
00215   h_vertex_sumTrks = bei->book1D("h_vertex_sumTrks", "Event Vertex, sum of track pt"      , 100,  0.0, 100.0 );
00216   h_vertex_d0    = bei->book1D("h_vertex_d0"   , "Event Vertex d0"                        , 100,  -10.0,  10.0);
00217   h_jet_et       = bei->book1D("h_jet_et",       "Jet with highest E_{T} (from "+theCaloJetCollectionLabel.label()+");E_{T}(1^{st} jet) (GeV)",    20, 0., 200.0);
00218   h_jet2_et      = bei->book1D("h_jet2_et",      "Jet with 2^{nd} highest E_{T} (from "+theCaloJetCollectionLabel.label()+");E_{T}(2^{nd} jet) (GeV)",    20, 0., 200.0);
00219   h_jet_count    = bei->book1D("h_jet_count",    "Number of "+theCaloJetCollectionLabel.label()+" (E_{T} > 15 GeV);Number of Jets", 8, -0.5, 7.5);
00220   h_caloMet      = bei->book1D("h_caloMet",        "Calo Missing E_{T}; GeV"          , 20,  0.0 , 100);
00221   h_caloMet_phi  = bei->book1D("h_caloMet_phi",    "Calo Missing E_{T} #phi;#phi(MET)", 35, -3.5, 3.5 );
00222   h_pfMet        = bei->book1D("h_pfMet",        "Pf Missing E_{T}; GeV"          , 20,  0.0 , 100);
00223   h_pfMet_phi    = bei->book1D("h_pfMet_phi",    "Pf Missing E_{T} #phi;#phi(MET)", 35, -3.5, 3.5 );
00224   h_eMultiplicity = bei_->book1D("NElectrons","# of electrons per event",10,0.,10.);
00225   h_mMultiplicity = bei_->book1D("NMuons","# of muons per event",10,0.,10.);
00226   h_ePt = bei_->book1D("ElePt","Pt of electrons",50,0.,100.);
00227   h_eEta = bei_->book1D("EleEta","Eta of electrons",100,-5.,5.);
00228   h_ePhi = bei_->book1D("ElePhi","Phi of electrons",100,-3.5,3.5);
00229   h_mPt_GMTM = bei_->book1D("MuonPt_GMTM","Pt of global+tracker muons",50,0.,100.);
00230   h_mEta_GMTM = bei_->book1D("MuonEta_GMTM","Eta of global+tracker muons",60,-3.,3.);
00231   h_mPhi_GMTM = bei_->book1D("MuonPhi_GMTM","Phi of global+tracker muons",70,-3.5,3.5);
00232   h_mPt_GMPT = bei_->book1D("MuonPt_GMPT","Pt of global prompt-tight muons",50,0.,100.);
00233   h_mEta_GMPT = bei_->book1D("MuonEta_GMPT","Eta of global prompt-tight muons",60,-3.,3.);
00234   h_mPhi_GMPT = bei_->book1D("MuonPhi_GMPT","Phi of global prompt-tight muons",70,-3.5,3.5);
00235   h_mPt_GM = bei_->book1D("MuonPt_GM","Pt of global muons",50,0.,100.);
00236   h_mEta_GM = bei_->book1D("MuonEta_GM","Eta of global muons",60,-3.,3.);
00237   h_mPhi_GM = bei_->book1D("MuonPhi_GM","Phi of global muons",70,-3.5,3.5);
00238   h_mPt_TM = bei_->book1D("MuonPt_TM","Pt of tracker muons",50,0.,100.);
00239   h_mEta_TM = bei_->book1D("MuonEta_TM","Eta of tracker muons",60,-3.,3.);
00240   h_mPhi_TM = bei_->book1D("MuonPhi_TM","Phi of tracker muons",70,-3.5,3.5);
00241   h_mPt_STAM = bei_->book1D("MuonPt_STAM","Pt of STA muons",50,0.,100.);
00242   h_mEta_STAM = bei_->book1D("MuonEta_STAM","Eta of STA muons",60,-3.,3.);
00243   h_mPhi_STAM = bei_->book1D("MuonPhi_STAM","Phi of STA muons",70,-3.5,3.5);
00244   h_eCombIso = bei_->book1D("EleCombIso","CombIso of electrons",100,0.,10.);
00245   h_mCombIso = bei_->book1D("MuonCombIso","CombIso of muons",100,0.,10.);
00246   h_dimumass_GMGM = bei->book1D("DimuMass_GMGM","Invariant mass of GMGM pairs",100,0.,200.);
00247   h_dimumass_GMTM = bei->book1D("DimuMass_GMTM","Invariant mass of GMTM pairs",100,0.,200.);
00248   h_dimumass_TMTM = bei->book1D("DimuMass_TMTM","Invariant mass of TMTM pairs",100,0.,200.);
00249   h_dielemass = bei->book1D("DieleMass","Invariant mass of EE pairs",100,0.,200.);
00250   h_lepcounts = bei->book1D("LeptonCounts","LeptonCounts for multi lepton events",49,0.,49.);
00251   
00252   bei->cd();  
00253 
00254 }
00255 //
00256 //  -- Analyze 
00257 //
00258 void HiggsDQM::analyze(const edm::Event& e, const edm::EventSetup& eSetup){
00259   //cout<<"[HiggsDQM::analyze()] "<<endl;
00260 
00261 //-------------------------------
00262 //--- Trigger Info
00263 //-------------------------------
00264   // short-circuit if hlt problems
00265   if( ! isValidHltConfig_ ) return;
00266   // Did it pass certain HLT path?
00267   Handle<TriggerResults> HLTresults;
00268   e.getByLabel(theTriggerResultsCollection, HLTresults); 
00269   if ( !HLTresults.isValid() ) return;
00270   //unsigned int triggerIndex_elec = hltConfig.triggerIndex(theElecTriggerPathToPass);
00271   //unsigned int triggerIndex_muon = hltConfig.triggerIndex(theMuonTriggerPathToPass);
00272   bool passed_electron_HLT = true;
00273   bool passed_muon_HLT     = true;
00274   //if (triggerIndex_elec < HLTresults->size()) passed_electron_HLT = HLTresults->accept(triggerIndex_elec);
00275   //if (triggerIndex_muon < HLTresults->size()) passed_muon_HLT     = HLTresults->accept(triggerIndex_muon);
00276   //if ( !(passed_electron_HLT || passed_muon_HLT) ) return;
00277 
00278 //-------------------------------
00279 //--- Vertex Info
00280 //-------------------------------
00281   Handle<VertexCollection> vertexHandle;
00282   e.getByLabel("offlinePrimaryVertices", vertexHandle);
00283   if ( vertexHandle.isValid() ){
00284     VertexCollection vertexCollection = *(vertexHandle.product());
00285     int vertex_number     = vertexCollection.size();
00286     VertexCollection::const_iterator v = vertexCollection.begin();
00287     double vertex_chi2    = v->normalizedChi2(); //v->chi2();
00288     double vertex_d0      = sqrt(v->x()*v->x()+v->y()*v->y());
00289     //double vertex_ndof    = v->ndof();cout << "ndof="<<vertex_ndof<<endl;
00290     double vertex_numTrks = v->tracksSize();
00291     double vertex_sumTrks = 0.0;
00292     for (Vertex::trackRef_iterator vertex_curTrack = v->tracks_begin(); vertex_curTrack!=v->tracks_end(); vertex_curTrack++) {
00293       vertex_sumTrks += (*vertex_curTrack)->pt();
00294     }
00295     h_vertex_number->Fill(vertex_number);
00296     h_vertex_chi2->Fill(vertex_chi2);
00297     h_vertex_d0  ->Fill(vertex_d0);
00298     h_vertex_numTrks->Fill(vertex_numTrks);
00299     h_vertex_sumTrks->Fill(vertex_sumTrks);
00300   }
00301   
00302 //-------------------------------
00303 //--- Electrons
00304 //-------------------------------
00305   float nEle=0;
00306   Handle<GsfElectronCollection> electronCollection;
00307   e.getByLabel(theElectronCollectionLabel, electronCollection);
00308   if ( electronCollection.isValid() ){
00309     int posEle=0,negEle=0;
00310     // If it passed electron HLT and the collection was found, find electrons near Z mass
00311     if( passed_electron_HLT ) {
00312       for (reco::GsfElectronCollection::const_iterator recoElectron=electronCollection->begin(); recoElectron!=electronCollection->end(); recoElectron++){
00313 //      cout << "Electron with pt= " <<  recoElectron->pt() << " and eta" << recoElectron->eta() << " p=" <<  recoElectron->p() << endl;
00314         h_ePt->Fill(recoElectron->pt());
00315         h_eEta->Fill(recoElectron->eta());
00316         h_ePhi->Fill(recoElectron->phi());
00317         if(recoElectron->charge()==1){
00318           posEle++;
00319         }else if(recoElectron->charge()==-1){
00320           negEle++;
00321         }
00322       // Require electron to pass some basic cuts
00323       //if ( recoElectron->et() < 20 || fabs(recoElectron->eta())>2.5 ) continue;
00324       // Tighter electron cuts
00325       //if ( recoElectron->deltaPhiSuperClusterTrackAtVtx() > 0.58 || 
00326       //     recoElectron->deltaEtaSuperClusterTrackAtVtx() > 0.01 || 
00327       //     recoElectron->sigmaIetaIeta() > 0.027 ) continue;
00328       } // end of loop over electrons
00329     } // end if passed HLT
00330     nEle = posEle+negEle; if(nEle>9.) nEle=9.;
00331     h_eMultiplicity->Fill(nEle);  
00332 
00333     // Z->ee:
00334     unsigned int eleCollectionSize = electronCollection->size();
00335     for(unsigned int i=0; i<eleCollectionSize; i++) {
00336       const GsfElectron& ele = electronCollection->at(i);
00337       double pt = ele.pt();
00338       if(pt>ptThrMu1_){
00339         for(unsigned int j=i+1; j<eleCollectionSize; j++) {
00340           const GsfElectron& ele2 = electronCollection->at(j);
00341           double pt2 = ele2.pt();
00342           if(pt2>ptThrMu2_){
00343             const math::XYZTLorentzVector ZRecoEE (ele.px()+ele2.px(), ele.py()+ele2.py() , ele.pz()+ele2.pz(), ele.p()+ele2.p());
00344             h_dielemass->Fill(ZRecoEE.mass());
00345           }
00346         }
00347       }
00348     }
00349   }
00350   
00351   
00352   
00353 //-------------------------------
00354 //--- Muons
00355 //-------------------------------
00356   float nMu=0;
00357   Handle<MuonCollection> muonCollection;
00358   e.getByLabel(theMuonCollectionLabel,muonCollection);
00359   if ( muonCollection.isValid() ){
00360     // Find the highest pt muons
00361     int posMu=0,negMu=0;
00362     TLorentzVector m1, m2;
00363     if( passed_muon_HLT ) {
00364       for (reco::MuonCollection::const_iterator recoMuon=muonCollection->begin(); recoMuon!=muonCollection->end(); recoMuon++){
00365         //cout << "Muon with pt= " <<  muIter->pt() << " and eta" << muIter->eta() << " p=" <<  muIter->p() << endl;
00366         if(recoMuon->isGlobalMuon()&&recoMuon->isTrackerMuon()){
00367           h_mPt_GMTM->Fill(recoMuon->pt());
00368           h_mEta_GMTM->Fill(recoMuon->eta());
00369           h_mPhi_GMTM->Fill(recoMuon->phi());
00370         }else if(recoMuon->isGlobalMuon()&&(muon::isGoodMuon( (*recoMuon),muon::GlobalMuonPromptTight))){
00371           h_mPt_GMPT->Fill(recoMuon->pt());
00372           h_mEta_GMPT->Fill(recoMuon->eta());
00373           h_mPhi_GMPT->Fill(recoMuon->phi());
00374         }else if(recoMuon->isGlobalMuon()){
00375           h_mPt_GM->Fill(recoMuon->pt());
00376           h_mEta_GM->Fill(recoMuon->eta());
00377           h_mPhi_GM->Fill(recoMuon->phi());
00378         }else if(recoMuon->isTrackerMuon()&&(muon::segmentCompatibility((*recoMuon),reco::Muon::SegmentAndTrackArbitration))){
00379           h_mPt_TM->Fill(recoMuon->pt());
00380           h_mEta_TM->Fill(recoMuon->eta());
00381           h_mPhi_TM->Fill(recoMuon->phi());
00382         }else if(recoMuon->isStandAloneMuon()){
00383           h_mPt_STAM->Fill(recoMuon->pt());
00384           h_mEta_STAM->Fill(recoMuon->eta());
00385           h_mPhi_STAM->Fill(recoMuon->phi());
00386         }
00387         if ( recoMuon->charge()==1 ){
00388           posMu++;
00389         }else if ( recoMuon->charge()==-1 ){
00390           negMu++;
00391         }
00392       }
00393       nMu = posMu+negMu; if(nMu>9.) nMu=9.;
00394       h_mMultiplicity->Fill(nMu);
00395     }
00396 
00397     // Z->mumu:
00398     unsigned int muonCollectionSize = muonCollection->size();
00399     for(unsigned int i=0; i<muonCollectionSize; i++) {
00400       const Muon& mu = muonCollection->at(i);
00401       //if (!mu.isGlobalMuon()) continue;
00402       double pt = mu.pt();
00403       if(pt>ptThrMu1_){
00404         for(unsigned int j=i+1; j<muonCollectionSize; j++) {
00405           const Muon& mu2 = muonCollection->at(j);
00406           double pt2 = mu2.pt();
00407           if(pt2>ptThrMu2_){
00408             // Glb + Glb  
00409             if(mu.isGlobalMuon() && mu2.isGlobalMuon()){
00410               const math::XYZTLorentzVector ZRecoGMGM (mu.px()+mu2.px(), mu.py()+mu2.py() , mu.pz()+mu2.pz(), mu.p()+mu2.p());
00411               h_dimumass_GMGM->Fill(ZRecoGMGM.mass());
00412             }
00413             // Glb + TM 
00414             else if(mu.isGlobalMuon() && mu2.isTrackerMuon()){
00415               const math::XYZTLorentzVector ZRecoGMTM (mu.px()+mu2.px(), mu.py()+mu2.py() , mu.pz()+mu2.pz(), mu.p()+mu2.p());
00416               h_dimumass_GMTM->Fill(ZRecoGMTM.mass());
00417             }
00418             // TM + TM 
00419             else if(mu.isTrackerMuon() && mu2.isTrackerMuon()){
00420               const math::XYZTLorentzVector ZRecoTMTM (mu.px()+mu2.px(), mu.py()+mu2.py() , mu.pz()+mu2.pz(), mu.p()+mu2.p());
00421               h_dimumass_TMTM->Fill(ZRecoTMTM.mass());
00422             }
00423           }
00424         }
00425       }
00426     }
00427   }
00428   
00429 //-------------------------------
00430 //--- Jets
00431 //-------------------------------
00432   Handle<CaloJetCollection> caloJetCollection;
00433   e.getByLabel (theCaloJetCollectionLabel,caloJetCollection);
00434   if ( caloJetCollection.isValid() ){
00435     float jet_et    = -8.0;
00436     //    float jet_eta   = -8.0; // UNUSED
00437     //    float jet_phi   = -8.0; // UNUSED
00438     int   jet_count = 0;
00439     float jet2_et   = -9.0;
00440     //    float jet2_eta  = -9.0; // UNUSED
00441     //    float jet2_phi  = -9.0; // UNUSED
00442     for (CaloJetCollection::const_iterator i_calojet = caloJetCollection->begin(); i_calojet != caloJetCollection->end(); i_calojet++) {
00443       float jet_current_et = i_calojet->et();
00444       // if it overlaps with electron, it is not a jet
00445       //if ( electron_et>0.0 && fabs(i_calojet->eta()-electron_eta ) < 0.2 && calcDeltaPhi(i_calojet->phi(), electron_phi ) < 0.2) continue;
00446       //if ( electron2_et>0.0&& fabs(i_calojet->eta()-electron2_eta) < 0.2 && calcDeltaPhi(i_calojet->phi(), electron2_phi) < 0.2) continue;
00447       // if it has too low Et, throw away
00448       if (jet_current_et < 15) continue;
00449       jet_count++;
00450       if (jet_current_et > jet_et) {
00451         jet2_et  = jet_et;  // 2nd highest jet get's et from current highest
00452         //        jet2_eta = jet_eta; // UNUSED
00453         //        jet2_phi = jet_phi; // UNUSED
00454         jet_et   = i_calojet->et(); // current highest jet gets et from the new highest
00455         //        jet_eta  = i_calojet->eta(); // UNUSED
00456         //        jet_phi  = i_calojet->phi(); // UNUSED
00457       } else if (jet_current_et > jet2_et) {
00458         jet2_et  = i_calojet->et();
00459         //        jet2_eta = i_calojet->eta(); // UNUSED
00460         //        jet2_phi = i_calojet->phi(); // UNUSED
00461       }
00462     }
00463     if (jet_et>0.0) {
00464       h_jet_et   ->Fill(jet_et);
00465       h_jet_count->Fill(jet_count);
00466     }
00467   }
00468   
00469 //-------------------------------
00470 //--- MET
00471 //-------------------------------
00472   Handle<CaloMETCollection> caloMETCollection;
00473   e.getByLabel(theCaloMETCollectionLabel, caloMETCollection);
00474   if ( caloMETCollection.isValid() ){
00475     float caloMet = caloMETCollection->begin()->et();
00476     float caloMet_phi = caloMETCollection->begin()->phi();
00477     h_caloMet        ->Fill(caloMet);
00478     h_caloMet_phi    ->Fill(caloMet_phi);
00479   }
00480   Handle<PFMETCollection> pfMETCollection;
00481   e.getByLabel(thePfMETCollectionLabel, pfMETCollection);
00482   if ( pfMETCollection.isValid() ){
00483     float pfMet = pfMETCollection->begin()->et();
00484     float pfMet_phi = pfMETCollection->begin()->phi();
00485     h_pfMet        ->Fill(pfMet);
00486     h_pfMet_phi    ->Fill(pfMet_phi);
00487   }
00488   
00489 //-------------------------------------
00490 //--- Events with more than 2 leptons:
00491 //-------------------------------------
00492   if(nMu+nEle > 2 && nMu+nEle < 10){
00493     if(nMu==0 && nEle==3) h_lepcounts->Fill(0);
00494     if(nMu==0 && nEle==4) h_lepcounts->Fill(1);
00495     if(nMu==0 && nEle==5) h_lepcounts->Fill(2);
00496     if(nMu==0 && nEle==6) h_lepcounts->Fill(3);
00497     if(nMu==0 && nEle==7) h_lepcounts->Fill(4);
00498     if(nMu==0 && nEle==8) h_lepcounts->Fill(5);
00499     if(nMu==0 && nEle==9) h_lepcounts->Fill(6);
00500     if(nMu==1 && nEle==2) h_lepcounts->Fill(7);
00501     if(nMu==1 && nEle==3) h_lepcounts->Fill(8);
00502     if(nMu==1 && nEle==4) h_lepcounts->Fill(9);
00503     if(nMu==1 && nEle==5) h_lepcounts->Fill(10);
00504     if(nMu==1 && nEle==6) h_lepcounts->Fill(11);
00505     if(nMu==1 && nEle==7) h_lepcounts->Fill(12);
00506     if(nMu==1 && nEle==8) h_lepcounts->Fill(13);
00507     if(nMu==2 && nEle==1) h_lepcounts->Fill(14);
00508     if(nMu==2 && nEle==2) h_lepcounts->Fill(15);
00509     if(nMu==2 && nEle==3) h_lepcounts->Fill(16);
00510     if(nMu==2 && nEle==4) h_lepcounts->Fill(17);
00511     if(nMu==2 && nEle==5) h_lepcounts->Fill(18);
00512     if(nMu==2 && nEle==6) h_lepcounts->Fill(19);
00513     if(nMu==2 && nEle==7) h_lepcounts->Fill(20);
00514     if(nMu==3 && nEle==0) h_lepcounts->Fill(21);
00515     if(nMu==3 && nEle==1) h_lepcounts->Fill(22);
00516     if(nMu==3 && nEle==2) h_lepcounts->Fill(23);
00517     if(nMu==3 && nEle==3) h_lepcounts->Fill(24);
00518     if(nMu==3 && nEle==4) h_lepcounts->Fill(25);
00519     if(nMu==3 && nEle==5) h_lepcounts->Fill(26);
00520     if(nMu==3 && nEle==6) h_lepcounts->Fill(27);
00521     if(nMu==4 && nEle==0) h_lepcounts->Fill(28);
00522     if(nMu==4 && nEle==1) h_lepcounts->Fill(29);
00523     if(nMu==4 && nEle==2) h_lepcounts->Fill(30);
00524     if(nMu==4 && nEle==3) h_lepcounts->Fill(31);
00525     if(nMu==4 && nEle==4) h_lepcounts->Fill(32);
00526     if(nMu==4 && nEle==5) h_lepcounts->Fill(33);
00527     if(nMu==5 && nEle==0) h_lepcounts->Fill(34);
00528     if(nMu==5 && nEle==1) h_lepcounts->Fill(35);
00529     if(nMu==5 && nEle==2) h_lepcounts->Fill(36);
00530     if(nMu==5 && nEle==3) h_lepcounts->Fill(37);
00531     if(nMu==5 && nEle==4) h_lepcounts->Fill(38);
00532     if(nMu==6 && nEle==0) h_lepcounts->Fill(39);
00533     if(nMu==6 && nEle==1) h_lepcounts->Fill(40);
00534     if(nMu==6 && nEle==2) h_lepcounts->Fill(41);
00535     if(nMu==6 && nEle==3) h_lepcounts->Fill(42);
00536     if(nMu==7 && nEle==0) h_lepcounts->Fill(43);
00537     if(nMu==7 && nEle==1) h_lepcounts->Fill(44);
00538     if(nMu==7 && nEle==2) h_lepcounts->Fill(45);
00539     if(nMu==8 && nEle==0) h_lepcounts->Fill(46);
00540     if(nMu==8 && nEle==1) h_lepcounts->Fill(47);
00541     if(nMu==9 && nEle==0) h_lepcounts->Fill(48);
00542     
00543   
00544   }
00545   if ((nMu+nEle) >= 10)
00546     LogDebug("HiggsDQM") <<"WARNING: "<<nMu+nEle<<" leptons in this event: run="<<e.id().run()<<", event="<<e.id().event()<< "\n"; 
00547 //     std::cout <<"WARNING: "<<nMu+nEle<<" leptons in this event: run="<<e.id().run()<<", event="<<e.id().event()<< "\n"; 
00548   
00549 /*  /// channel conditions
00550   nElectron=0;
00551   nMuon=0;
00552   if (decaychannel=="2e2mu") {
00553     if ( posEle>=1 && negEle>=1 ) {
00554       nElectron=posEle+negEle;
00555     }
00556     if ( posMu>=1 && negMu>=1 ) {
00557       nMuon=posMu+negMu;
00558     }
00559   }
00560   else if (decaychannel=="4e") {
00561     if ( posEle>=2 && negEle>=2 ) {
00562       nElectron=posEle+negEle;
00563     }
00564   }
00565   else if (decaychannel=="4mu") {
00566     if ( posMu>=2 && negMu>=2 ) {
00567       nMuon=posMu+negMu;
00568     }  
00569   }
00570        
00571 
00572   // Pairs of EE MuMu 
00573   int nZEE=0,nZMuMu=0;
00574   if (decaychannel=="2e2mu"){
00575     Handle<CompositeCandidateCollection> zEECandidates;
00576     e.getByLabel(zToEETag_.label(), zEECandidates);    
00577     for ( CompositeCandidateCollection::const_iterator zIter=zEECandidates->begin(); zIter!= zEECandidates->end(); ++zIter ) {
00578 //      cout << "Zee mass= " << double(zIter->p4().mass()) << endl;
00579       if ( double(zIter->p4().mass())> 12. ){ 
00580         nZEE++;
00581       }
00582     }
00583     Handle<CompositeCandidateCollection> zMuMuCandidates;
00584     e.getByLabel(zToMuMuTag_.label(), zMuMuCandidates);
00585     for ( CompositeCandidateCollection::const_iterator zIter=zMuMuCandidates->begin(); zIter!= zMuMuCandidates->end(); ++zIter ) {
00586 //      cout << "Zmumu mass= " << double(zIter->p4().mass()) << endl;
00587       if ( zIter->p4().mass()> 12. ){
00588         nZMuMu++;
00589       }
00590     }    
00591   }
00592   
00593   // exclusive couples ZEE and ZMUMU
00594   if (decaychannel=="4e" ){
00595     Handle<CompositeCandidateCollection> higgsCandidates;
00596     e.getByLabel(hTozzTo4leptonsTag_.label(), higgsCandidates);
00597     for ( CompositeCandidateCollection::const_iterator hIter=higgsCandidates->begin(); hIter!= higgsCandidates->end(); ++hIter ) {
00598       if (nZEE<2) nZEE=0;
00599       for (size_t ndau=0; ndau<hIter->numberOfDaughters();ndau++){
00600         if ( hIter->daughter(ndau)->p4().mass()> 12.){  
00601           nZEE++;
00602         }
00603       }
00604     }
00605   }
00606   //
00607   if (decaychannel=="4mu" ){
00608     Handle<CompositeCandidateCollection> higgsCandidates;
00609     e.getByLabel(hTozzTo4leptonsTag_.label(), higgsCandidates);
00610     for ( CompositeCandidateCollection::const_iterator hIter=higgsCandidates->begin(); hIter!= higgsCandidates->end(); ++hIter ) {
00611       if (nZMuMu<2) nZMuMu=0;
00612       for (size_t ndau=0; ndau<hIter->numberOfDaughters();ndau++){
00613         if ( hIter->daughter(ndau)->p4().mass()> 12. ){  
00614           nZMuMu++;
00615         }
00616       }
00617     }
00618   }
00619 
00620   // 4 lepton combinations
00621   Handle<CompositeCandidateCollection> higgsCandidates;
00622   e.getByLabel(hTozzTo4leptonsTag_.label(), higgsCandidates);
00623    int nHiggs=0;
00624   for ( CompositeCandidateCollection::const_iterator hIter=higgsCandidates->begin(); hIter!= higgsCandidates->end(); ++hIter ) {
00625     if ( hIter->p4().mass()> 100. ){  
00626       nHiggs++;
00627     }
00628   }    
00629 */
00630 //  cout<<"Z and Higgs candidates: "<<nZEE<<" "<<nZMuMu<<" "<<nHiggs<<endl;
00631 
00632   //cout<<"[leaving HiggsDQM::analyze()] "<<endl;
00633 
00634 }
00635 //
00636 // -- End Luminosity Block
00637 //
00638 void HiggsDQM::endLuminosityBlock(edm::LuminosityBlock const& lumiSeg, edm::EventSetup const& eSetup) {
00639 //  cout<<"Entering HiggsDQM::endLuminosityBlock: "<<endl;
00640 
00641   edm::LogInfo ("HiggsDQM") <<"[HiggsDQM]: End of LS transition, performing the DQM client operation";
00642 
00643   nLumiSecs_++;
00644   //cout << "nLumiSecs_: "<< nLumiSecs_ << endl;
00645   
00646   edm::LogInfo("HiggsDQM") << "====================================================== " << endl << " ===> Iteration # " << nLumiSecs_ << " " << lumiSeg.luminosityBlock() << endl  << "====================================================== " << endl;
00647 
00648 //  cout<<"...leaving HiggsDQM::endLuminosityBlock. "<<endl;
00649 }
00650 //
00651 // -- End Run
00652 //
00653 void HiggsDQM::endRun(edm::Run const& run, edm::EventSetup const& eSetup){
00654 //  cout<<"Entering HiggsDQM::endRun: "<<endl;
00655 
00656   //edm::LogVerbatim ("HiggsDQM") <<"[HiggsDQM]: End of Run, saving  DQM output ";
00657   //int iRun = run.run();
00658   
00659 //  cout<<"...leaving HiggsDQM::endRun. "<<endl;
00660 }
00661 
00662 //
00663 // -- End Job
00664 //
00665 void HiggsDQM::endJob(){
00666 //  cout<<"In HiggsDQM::endJob "<<endl;
00667   edm::LogInfo("HiggsDQM") <<"[HiggsDQM]: endjob called!";
00668 
00669 }
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/DQM/Physics/src/HiggsDQM.cc
