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