Reference Manual

Go to the documentation of this file.
00001 /*
00002  *  See header file for a description of this class.
00003  *
00004  *  $Date: 2012/06/19 10:13:04 $
00005  *  $Revision: 1.20 $
00006  *  \author Valentina Gori, University of Firenze
00007  */
00008 
00009 #include "DQM/Physics/src/EwkDQM.h"
00010 #include "DataFormats/Candidate/interface/Candidate.h"
00011 #include "FWCore/Framework/interface/Event.h"
00012 #include "FWCore/Framework/interface/EventSetup.h"
00013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00014 
00015 #include "DQMServices/Core/interface/DQMStore.h"
00016 #include "DQMServices/Core/interface/MonitorElement.h"
00017 #include "FWCore/ServiceRegistry/interface/Service.h"
00018 
00019 #include "DataFormats/Common/interface/Handle.h"
00020 
00021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00022 #include "FWCore/Common/interface/TriggerNames.h"
00023 
00024 // Physics Objects
00025 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
00026 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
00027 #include "DataFormats/MuonReco/interface/Muon.h"
00028 #include "DataFormats/MuonReco/interface/MuonFwd.h"
00029 #include "DataFormats/JetReco/interface/Jet.h"
00030 #include "DataFormats/METReco/interface/MET.h"
00031 #include "DataFormats/TrackReco/interface/Track.h"
00032 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00033 #include "DataFormats/VertexReco/interface/Vertex.h"
00034 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00035 
00036 
00037 #include "DataFormats/Math/interface/LorentzVector.h"
00038 #include "TLorentzVector.h"
00039 
00040 #include <vector>
00041 
00042 #include <string>
00043 #include <cmath>
00044 using namespace std;
00045 using namespace edm;
00046 using namespace reco;
00047 
00048 
00049 
00050 // EwkDQM::EwkDQM(const ParameterSet& parameters) {
00051 EwkDQM::EwkDQM(const ParameterSet& parameters) {
00052   eJetMin_     = parameters.getUntrackedParameter<double>("EJetMin", 999999.);
00053 
00054   // riguardare questa sintassi 
00055   // Get parameters from configuration file
00056   theElecTriggerPathToPass_   = parameters.getParameter<string>("elecTriggerPathToPass");
00057   theMuonTriggerPathToPass_    = parameters.getParameter<string>("muonTriggerPathToPass");
00058   //eleTrigPathNames_ = parameters.getUntrackedParameter< std::vector<std::string> >("eleTrigPathNames");
00059   //muTrigPathNames_ = parameters.getUntrackedParameter< std::vector<std::string> >("muTrigPathNames");
00060   theTriggerResultsCollection_ = parameters.getParameter<InputTag>("triggerResultsCollection");
00061   theMuonCollectionLabel_      = parameters.getParameter<InputTag>("muonCollection");
00062   theElectronCollectionLabel_  = parameters.getParameter<InputTag>("electronCollection");
00063   //  theCaloJetCollectionLabel_   = parameters.getParameter<InputTag>("caloJetCollection");
00064   thePFJetCollectionLabel_   = parameters.getParameter<InputTag>("PFJetCollection");
00065   theCaloMETCollectionLabel_   = parameters.getParameter<InputTag>("caloMETCollection");
00066   
00067   // just to initialize
00068   isValidHltConfig_ = false;
00069 
00070   // coverity says.. (cos'e` questo?? che variabili sono???)
00071   h_vertex_number = 0;
00072   h_vertex_chi2 = 0;
00073   h_vertex_numTrks = 0;
00074   h_vertex_sumTrks = 0;
00075   h_vertex_d0 = 0;
00076 
00077   h_jet_count = 0;
00078   h_jet_et = 0;
00079   h_jet_pt = 0;   // prova
00080   h_jet_eta = 0;  // aggiunto il 23 maggio 
00081   h_jet_phi = 0;
00082   h_jet2_et = 0;
00083   //h_jet2_pt = 0;
00084   h_jet2_eta = 0;
00085   h_jet2_phi = 0;
00086 
00087   h_e1_et = 0;
00088   h_e2_et = 0;
00089   h_e1_eta = 0;
00090   h_e2_eta = 0;
00091   h_e1_phi = 0;
00092   h_e2_phi = 0;
00093 
00094   h_m1_pt = 0;
00095   h_m2_pt = 0;
00096   h_m1_eta = 0;
00097   h_m2_eta = 0;
00098   h_m1_phi = 0;
00099   h_m2_phi = 0;
00100 
00101   //h_t1_et = 0;
00102   //h_t1_eta = 0;
00103   //h_t1_phi = 0;
00104 
00105   h_met = 0;
00106   h_met_phi = 0;
00107 
00108   h_e_invWMass = 0;
00109   h_m_invWMass = 0;
00110   h_mumu_invMass = 0;
00111   h_ee_invMass = 0;
00112   
00113   theDbe = Service<DQMStore>().operator->();
00114 
00115 }
00116 
00117 EwkDQM::~EwkDQM() { 
00118 }
00119 
00120 
00121 void EwkDQM::beginJob() {
00122 
00123   char chtitle[256] = "";
00124 
00125   logTraceName = "EwkAnalyzer";
00126 
00127   LogTrace(logTraceName)<<"Parameters initialization";
00128   theDbe->setCurrentFolder("Physics/EwkDQM");  // Use folder with name of PAG
00129 
00130   const float pi = 4*atan(1);
00131 
00132   // Keep the number of plots and number of bins to a minimum!
00133   h_vertex_number = theDbe->book1D("vertex_number", "Number of event vertices in collection", 10,-0.5,   9.5 );
00134   h_vertex_chi2  = theDbe->book1D("vertex_chi2" , "Event Vertex #chi^{2}/n.d.o.f."          , 20, 0.0,   2.0 );
00135   h_vertex_numTrks = theDbe->book1D("vertex_numTrks", "Event Vertex, number of tracks"     , 20, -0.5,  59.5 );
00136   h_vertex_sumTrks = theDbe->book1D("vertex_sumTrks", "Event Vertex, sum of track pt"      , 20,  0.0, 100.0 );
00137   h_vertex_d0    = theDbe->book1D("vertex_d0"   , "Event Vertex d0"                        , 20,  0.0,   0.05);
00138   
00139   snprintf(chtitle, 255,  "Number of %s (E_{T} > 15 GeV);Number of Jets",thePFJetCollectionLabel_.label().data());
00140   h_jet_count    = theDbe->book1D("jet_count", chtitle, 8, -0.5, 7.5);
00141 
00142   snprintf(chtitle, 255, "Leading jet E_{T} (from %s);E_{T}(1^{st} jet) (GeV)",
00143            thePFJetCollectionLabel_.label().data());
00144   h_jet_et       = theDbe->book1D("jet_et", chtitle,    20, 0., 200.0);
00145 
00146   snprintf(chtitle, 255, "Leading jet p_{T} (from %s);p_{T}(1^{st} jet) (GeV/c)", thePFJetCollectionLabel_.label().data());
00147   h_jet_pt       = theDbe->book1D("jet_pt", chtitle,  20, 0., 200.0); 
00148 
00149   snprintf(chtitle, 255,  "Leading jet #eta (from %s); #eta (1^{st} jet)",thePFJetCollectionLabel_.label().data());
00150   h_jet_eta      = theDbe->book1D("jet_eta", chtitle,  20, -10., 10.0);
00151   snprintf(chtitle, 255, "Leading jet #phi (from %s); #phi(1^{st} jet)", thePFJetCollectionLabel_.label().data());
00152   h_jet_phi      = theDbe->book1D("jet_phi", chtitle,  22, -1.1*pi, 1.1*pi); 
00153 
00154   snprintf(chtitle, 255, "2^{nd} leading jet E_{T} (from %s);E_{T}(2^{nd} jet) (GeV)",thePFJetCollectionLabel_.label().data());
00155   h_jet2_et      = theDbe->book1D("jet2_et", chtitle,  20, 0., 200.0);
00156   //snprintf(chtitle, 255, "2^{nd} leading jet p_{T} (from %s);p_{T}(2^{nd} jet) (GeV/c)", thePFJetCollectionLabel_.label().data());
00157   //h_jet2_pt       = theDbe->book1D("jet2_pt", chtitle,  20, 0., 200.0); 
00158   snprintf(chtitle, 255,  "2^{nd} leading jet #eta (from %s); #eta (2^{nd} jet)",thePFJetCollectionLabel_.label().data());
00159   h_jet2_eta      = theDbe->book1D("jet2_eta", chtitle,  20, -10., 10.0);
00160   snprintf(chtitle, 255, "2^{nd} leading jet #phi (from %s); #phi(2^{nd} jet)", thePFJetCollectionLabel_.label().data());
00161   h_jet2_phi      = theDbe->book1D("jet2_phi", chtitle,  22, -1.1*pi, 1.1*pi); 
00162 
00163   h_e1_et        = theDbe->book1D("e1_et",  "E_{T} of Leading Electron;E_{T} (GeV)"        , 20,  0.0 , 100.0);
00164   h_e2_et        = theDbe->book1D("e2_et",  "E_{T} of Second Electron;E_{T} (GeV)"         , 20,  0.0 , 100.0);
00165   h_e1_eta       = theDbe->book1D("e1_eta", "#eta of Leading Electron;#eta"                , 20, -4.0 , 4.0);
00166   h_e2_eta       = theDbe->book1D("e2_eta", "#eta of Second Electron;#eta"                 , 20, -4.0 , 4.0);
00167   h_e1_phi       = theDbe->book1D("e1_phi", "#phi of Leading Electron;#phi"                , 22, -1.1*pi, 1.1*pi );
00168   h_e2_phi       = theDbe->book1D("e2_phi", "#phi of Second Electron;#phi"                 , 22, -1.1*pi, 1.1*pi );
00169   h_m1_pt        = theDbe->book1D("m1_pt",  "p_{T} of Leading Muon;p_{T}(1^{st} #mu) (GeV)", 20,  0.0 , 100.0);
00170   h_m2_pt        = theDbe->book1D("m2_pt",  "p_{T} of Second Muon;p_{T}(2^{nd} #mu) (GeV)" , 20,  0.0 , 100.0);
00171   h_m1_eta       = theDbe->book1D("m1_eta", "#eta of Leading Muon;#eta(1^{st} #mu)"        , 20, -4.0 , 4.0);
00172   h_m2_eta       = theDbe->book1D("m2_eta", "#eta of Second Muon;#eta(2^{nd} #mu)"         , 20, -4.0 , 4.0);
00173   h_m1_phi       = theDbe->book1D("m1_phi", "#phi of Leading Muon;#phi(1^{st} #mu)"        , 20, (-1.-1./10.)*pi, (1.+1./10.)*pi);
00174   h_m2_phi       = theDbe->book1D("m2_phi", "#phi of Second Muon;#phi(2^{nd} #mu)"         , 20, (-1.-1./10.)*pi, (1.+1./10.)*pi);
00175 //  h_t1_et          = theDbe->book1D("t1_et",           "E_{T} of Leading Tau;E_{T} (GeV)" , 20, 0.0 , 100.0);
00176 //  h_t1_eta         = theDbe->book1D("t1_eta",          "#eta of Leading Tau;#eta"               , 20, -4.0, 4.0);
00177 //  h_t1_phi         = theDbe->book1D("t1_phi",          "#phi of Leading Tau;#phi"               , 20, -4.0, 4.0);
00178   snprintf(chtitle, 255, "Missing E_{T} (%s); GeV", theCaloMETCollectionLabel_.label().data());
00179   h_met          = theDbe->book1D("met",  chtitle, 20,  0.0 , 100);
00180   h_met_phi      = theDbe->book1D("met_phi",    "Missing E_{T} #phi;#phi(MET)"             , 22, (-1.-1./10.)*pi, (1.+1./10.)*pi );
00181 
00182   h_e_invWMass   = theDbe->book1D("we_invWMass", "W-> e #nu Transverse Mass;M_{T} (GeV)"    , 20,  0.0, 140.0); 
00183   h_m_invWMass   = theDbe->book1D("wm_invWMass", "W-> #mu #nu Transverse Mass;M_{T} (GeV)"  , 20,  0.0, 140.0); 
00184 
00185   h_mumu_invMass = theDbe->book1D("z_mm_invMass", "#mu#mu Invariant Mass;InvMass (GeV)"    , 20, 40.0, 140.0 );
00186   h_ee_invMass   = theDbe->book1D("z_ee_invMass",   "ee Invariant Mass;InvMass (Gev)"        , 20, 40.0, 140.0 );
00187 }
00188 
00189 
00190 
00194 void EwkDQM::beginRun( const edm::Run& theRun, const edm::EventSetup& theSetup ) {
00195   
00196   // passed as parameter to HLTConfigProvider::init(), not yet used
00197   bool isConfigChanged = false;
00198   
00199   // isValidHltConfig_ used to short-circuit analyze() in case of problems
00200   const std::string hltProcessName( theTriggerResultsCollection_.process() );
00201   isValidHltConfig_ = hltConfigProvider_.init( theRun, theSetup, hltProcessName, isConfigChanged );
00202 
00203 }
00204 
00205 
00206 void EwkDQM::analyze(const Event& iEvent, const EventSetup& iSetup) {
00207 
00208   // short-circuit if hlt problems
00209   if( ! isValidHltConfig_ ) return;
00210 
00211   // non mi e` chiaro come faccia a "ciclare" sugli eventi
00212   LogTrace(logTraceName)<<"Analysis of event # ";
00213   // Did it pass certain HLT path?
00214   Handle<TriggerResults> HLTresults;
00215   iEvent.getByLabel(theTriggerResultsCollection_, HLTresults); 
00216   if ( !HLTresults.isValid() ) return;
00217 
00218   const edm::TriggerNames & trigNames = iEvent.triggerNames(*HLTresults);
00219 
00220   // a temporary, until we have a list of triggers of interest
00221   std::vector<std::string> eleTrigPathNames;
00222   std::vector<std::string> muTrigPathNames;
00223   eleTrigPathNames.push_back(theElecTriggerPathToPass_);
00224   muTrigPathNames.push_back(theMuonTriggerPathToPass_);
00225   // end of temporary
00226 
00227   bool passed_electron_HLT = false;
00228   bool passed_muon_HLT     = false;
00229   for (unsigned int i=0; i<HLTresults->size(); i++) {
00230     const std::string trigName = trigNames.triggerName(i);
00231     // check if triggerName matches electronPath
00232     for(unsigned int index=0; index<eleTrigPathNames.size() && !passed_electron_HLT; index++) {
00233       size_t trigPath = trigName.find(eleTrigPathNames[index]); // 0 if found, pos if not
00234       if (trigPath==0) {
00235         passed_electron_HLT = HLTresults->accept(i);
00236       }
00237     }
00238     // check if triggerName matches muonPath
00239     for(unsigned int index=0; index<muTrigPathNames.size() && !passed_muon_HLT; index++) {
00240       size_t trigPath = trigName.find(muTrigPathNames[index]); // 0 if found, pos if not
00241       if (trigPath==0) {
00242         passed_muon_HLT = HLTresults->accept(i);
00243       }
00244     }
00245   }
00246   
00247   // we are interested in events with a valid electron or muon
00248   if ( !(passed_electron_HLT || passed_muon_HLT) ) return;
00249 
00251   //Vertex information
00252   Handle<VertexCollection> vertexHandle;
00253   iEvent.getByLabel("offlinePrimaryVertices", vertexHandle);
00254   if ( !vertexHandle.isValid() ) return;
00255   VertexCollection vertexCollection = *(vertexHandle.product());
00256   int vertex_number     = vertexCollection.size();
00257   VertexCollection::const_iterator v = vertexCollection.begin();
00258   double vertex_chi2    = v->normalizedChi2(); //v->chi2();
00259   double vertex_d0      = sqrt(v->x()*v->x()+v->y()*v->y());
00260   //std::cout << "vertex_d0=" << vertex_d0 << "\n";
00261   //double vertex_ndof    = v->ndof();cout << "ndof="<<vertex_ndof<<endl;
00262   double vertex_numTrks = v->tracksSize();
00263   double vertex_sumTrks = 0.0;
00264   for (Vertex::trackRef_iterator vertex_curTrack = v->tracks_begin(); vertex_curTrack!=v->tracks_end(); vertex_curTrack++) {
00265     vertex_sumTrks += (*vertex_curTrack)->pt();
00266   }
00267 
00269   //Missing ET
00270   Handle< View<MET> > caloMETCollection;
00271   iEvent.getByLabel(theCaloMETCollectionLabel_, caloMETCollection);
00272   if ( !caloMETCollection.isValid() ) return;
00273   float missing_et = caloMETCollection->begin()->et();
00274   float met_phi    = caloMETCollection->begin()->phi();
00275 
00276 
00278   // grab "gaussian sum fitting" electrons
00279   Handle<GsfElectronCollection> electronCollection;
00280   iEvent.getByLabel(theElectronCollectionLabel_, electronCollection);
00281   if ( !electronCollection.isValid() ) return;
00282 
00283   // Find the highest and 2nd highest electron
00284   float electron_et   = -8.0;
00285   float electron_eta  = -8.0;
00286   float electron_phi  = -8.0;
00287   float electron2_et  = -9.0;
00288   float electron2_eta = -9.0;
00289   float electron2_phi = -9.0;
00290   float ee_invMass    = -9.0;
00291   TLorentzVector e1, e2;
00292 
00293   // If it passed electron HLT and the collection was found, find electrons near Z mass
00294   if( passed_electron_HLT ) {
00295 
00296     for (reco::GsfElectronCollection::const_iterator recoElectron=electronCollection->begin(); recoElectron!=electronCollection->end(); recoElectron++){
00297 
00298       // Require electron to pass some basic cuts
00299       if ( recoElectron->et() < 20 || fabs(recoElectron->eta())>2.5 ) continue;
00300 
00301       // Tighter electron cuts
00302       if ( recoElectron->deltaPhiSuperClusterTrackAtVtx() > 0.58 || 
00303            recoElectron->deltaEtaSuperClusterTrackAtVtx() > 0.01 || 
00304            recoElectron->sigmaIetaIeta() > 0.027 ) continue;
00305 
00306       if (recoElectron->et() > electron_et){
00307         electron2_et  = electron_et;  // 2nd highest gets values from current highest
00308         electron2_eta = electron_eta;
00309         electron2_phi = electron_phi;
00310         electron_et   = recoElectron->et();  // 1st highest gets values from new highest
00311         electron_eta  = recoElectron->eta();
00312         electron_phi  = recoElectron->phi();
00313         e1 = TLorentzVector(recoElectron->momentum().x(),recoElectron->momentum().y(),recoElectron->momentum().z(),recoElectron->p());
00314       } else if (recoElectron->et() > electron2_et) {
00315         electron2_et  = recoElectron->et();
00316         electron2_eta = recoElectron->eta();
00317         electron2_phi = recoElectron->phi();
00318         e2 = TLorentzVector(recoElectron->momentum().x(),recoElectron->momentum().y(),recoElectron->momentum().z(),recoElectron->p());
00319       }
00320     } // end of loop over electrons
00321     if (electron2_et>0.0) {
00322       TLorentzVector pair=e1+e2;
00323       ee_invMass = pair.M();
00324     }
00325   } // end of "are electrons valid"
00327 
00328 
00329 
00331   // Take the STA muon container
00332   Handle<MuonCollection> muonCollection;
00333   iEvent.getByLabel(theMuonCollectionLabel_,muonCollection);
00334   if ( !muonCollection.isValid() ) return;
00335 
00336   // Find the highest pt muons
00337   float mm_invMass = -9.0;
00338   float muon_pt   = -9.0;
00339   float muon_eta  = -9.0;
00340   float muon_phi  = -9.0;
00341   float muon2_pt  = -9.0;
00342   float muon2_eta = -9.0;
00343   float muon2_phi = -9.0;
00344   TLorentzVector m1, m2;
00345 
00346   if( passed_muon_HLT ) {
00347     for (reco::MuonCollection::const_iterator recoMuon=muonCollection->begin(); recoMuon!=muonCollection->end(); recoMuon++){
00348 
00349       // Require muon to pass some basic cuts
00350       if ( recoMuon->pt() < 20 || !recoMuon->isGlobalMuon() ) continue;
00351       // Some tighter muon cuts
00352       if ( recoMuon->globalTrack()->normalizedChi2() > 10 ) continue;
00353 
00354       if (recoMuon->pt() > muon_pt){
00355         muon2_pt  = muon_pt;  // 2nd highest gets values from current highest    
00356         muon2_eta = muon_eta;
00357         muon2_phi = muon_phi;
00358         muon_pt   = recoMuon->pt();  // 1st highest gets values from new highest
00359         muon_eta  = recoMuon->eta();
00360         muon_phi  = recoMuon->phi();
00361         m1 = TLorentzVector(recoMuon->momentum().x(),recoMuon->momentum().y(),recoMuon->momentum().z(),recoMuon->p());
00362       } else if (recoMuon->pt() > muon2_pt) {
00363         muon2_pt  = recoMuon->pt();
00364         muon2_eta = recoMuon->eta();
00365         muon2_phi = recoMuon->phi();
00366         m2 = TLorentzVector(recoMuon->momentum().x(),recoMuon->momentum().y(),recoMuon->momentum().z(),recoMuon->p());
00367       }
00368     }
00369   }
00370   if (muon2_pt>0.0) {
00371     TLorentzVector pair=m1+m2;
00372     mm_invMass = pair.M();
00373   }
00375   
00376   
00378   // Find the highest et jet                    
00379     
00380   //  Handle<CaloJetCollection> caloJetCollection;
00381   Handle<View<Jet> > PFJetCollection;
00382   //  iEvent.getByLabel (theCaloJetCollectionLabel,caloJetCollection);
00383   iEvent.getByLabel (thePFJetCollectionLabel_,PFJetCollection);
00384   //  if ( !caloJetCollection.isValid() ) return;
00385   if ( !PFJetCollection.isValid() ) return;
00386   
00387   unsigned int muonCollectionSize = muonCollection->size();
00388   //unsigned int jetCollectionSize = jetCollection->size();
00389   unsigned int PFJetCollectionSize = PFJetCollection->size();
00390   int jet_count = 0; 
00391   //int LEADJET=-1;  double max_pt=0;
00392   
00393   
00394   float jet_et    = -80.0;
00395   float jet_pt    = -80.0;  // prova
00396   float jet_eta   = -80.0; // now USED
00397   float jet_phi   = -80.0; // now USED
00398   float jet2_et   = -90.0;
00399   float jet2_eta  = -90.0; // now USED
00400   float jet2_phi  = -90.0; // now USED
00401   //  for (CaloJetCollection::const_iterator i_calojet = caloJetCollection->begin(); i_calojet != caloJetCollection->end(); i_calojet++) {
00402   //  for (PFJetCollection::const_iterator i_pfjet = PFJetCollection->begin(); i_pfjet != PFJetCollection->end(); i_pfjet++) {
00403 
00404 
00405 
00406 
00407    //  float jet_current_et = i_calojet->et();
00408    //  float jet_current_et = i_pfjet->et();            // e` identico a jet.et()
00409 
00410 
00411 
00412     //    jet_count++;
00413 
00414 
00415 
00416 
00417     // cleaning: va messo prima del riempimento dell'istogramma // This is in order to use PFJets
00418     for (unsigned int i=0; i<PFJetCollectionSize; i++) {
00419       const Jet& jet = PFJetCollection->at(i);
00420       // la classe "jet" viene definita qui!!!
00421       double minDistance=99999;
00422       for (unsigned int j=0; j<muonCollectionSize; j++) {
00423         const Muon& mu = muonCollection->at(j);
00424         double distance = sqrt( (mu.eta()-jet.eta())*(mu.eta()-jet.eta()) +(mu.phi()-jet.phi())*(mu.phi()-jet.phi()) );      
00425         if (minDistance>distance) minDistance=distance;
00426       }
00427       if (minDistance<0.3) continue; // 0.3 is the isolation cone around the muon 
00428       // se la distanza muone-cono del jet e` minore di 0.3, passo avanti e non conteggio il mio jet
00429        
00430       // If it overlaps with ELECTRON, it is not a jet
00431       if ( electron_et>0.0 && fabs(jet.eta()-electron_eta ) < 0.2 && calcDeltaPhi(jet.phi(), electron_phi ) < 0.2) continue;
00432       if ( electron2_et>0.0&& fabs(jet.eta()-electron2_eta) < 0.2 && calcDeltaPhi(jet.phi(), electron2_phi) < 0.2) continue;
00433       
00434       // provo a cambiare la parte degli elettroni in modo simmetrico alla parte per i muoni
00435 
00436       // ...
00437       // ...
00438 
00439       
00440       // if it has too low Et, throw away
00441       if (jet.et() < eJetMin_) continue;
00442       jet_count ++;
00443 
00444       // ovvero: incrementa jet_count se:
00445       //   - non c'e un muone entro 0.3 di distanza dal cono del jet;
00446       //   - se il jet non si sovrappone ad un elettrone;
00447       //   - se l'energia trasversa e` maggiore della soglia impostata (15?)
00448       
00449       //if(jet.et()>max_pt) { LEADJET=i; max_pt=jet.et();} 
00450       // se l'energia del jet e` maggiore di max_pt, diventa "i" l'indice del jet piu` energetico e max_pt la sua energia 
00451         
00452   
00453     // riguardare questo!!!
00454     // fino ad ora, jet_et era inizializzato a -8.0
00455     if (jet.et() > jet_et) {
00456       jet2_et  = jet_et;  // 2nd highest jet gets et from current highest // perche` prende l'energia del primo jet??
00457       jet2_eta = jet_eta; // now USED   
00458       jet2_phi = jet_phi; // now USED
00459       //      jet_et   = i_calojet->et(); // current highest jet gets et from the new highest
00460       jet_et   = jet.et(); // current highest jet gets et from the new highest
00461       // ah, ok! lo riaggiorna solo dopo!
00462       jet_pt = jet.pt();           // e` il pT del leading jet
00463       jet_eta  = jet.eta(); // now USED
00464       jet_phi  = jet.phi()*(Geom::pi()/180.); // now USED
00465     } 
00466         else if (jet.et() > jet2_et) {
00467           //      jet2_et  = i_calojet->et();
00468           jet2_et  = jet.et();
00469           //      jet2_eta = i_calojet->eta(); // UNUSED
00470           //      jet2_phi = i_calojet->phi(); // UNUSED
00471           jet2_eta = jet.eta(); // now USED
00472           jet2_phi = jet.phi(); // now USED
00473         }
00474     // questo elseif funziona
00475     
00476     }
00478 
00479 
00480 
00482   //                 Fill Histograms                                            //
00484 
00485   bool fill_e1  = false;
00486   bool fill_e2  = false;
00487   bool fill_m1  = false;
00488   bool fill_m2  = false;
00489   bool fill_met = false;
00490 
00491   // Was Z->ee found?
00492   if (ee_invMass>0.0) {
00493     h_ee_invMass ->Fill(ee_invMass);
00494     fill_e1 = true;
00495     fill_e2 = true;
00496   }
00497 
00498   // Was Z->mu mu found?
00499   if (mm_invMass > 0.0) {
00500     h_mumu_invMass->Fill(mm_invMass);
00501     fill_m1 = true;
00502     fill_m2 = true;    h_jet2_et  ->Fill(jet2_et);
00503   }
00504 
00505   // Was W->e nu found?
00506   if (electron_et>0.0&&missing_et>20.0) {
00507     float dphiW = fabs(met_phi-electron_phi);
00508     float W_mt_e = sqrt(2*missing_et*electron_et*(1-cos(dphiW)));
00509     h_e_invWMass ->Fill(W_mt_e);
00510     fill_e1  = true;
00511     fill_met = true;
00512   }
00513 
00514   // Was W->mu nu found?
00515   if (muon_pt>0.0&&missing_et>20.0) {
00516     float dphiW = fabs(met_phi-muon_phi);
00517     float W_mt_m = sqrt(2*missing_et*muon_pt*(1-cos(dphiW)));
00518     h_m_invWMass ->Fill(W_mt_m);
00519     fill_m1  = true;
00520     fill_met = true;
00521   }
00522 
00523 
00524 
00525   if (jet_et>-10.0) {
00526     h_jet_et   ->Fill(jet_et);
00527     h_jet_count->Fill(jet_count);
00528   }
00529 
00530   if (jet_pt>0.) {
00531     h_jet_pt   ->Fill(jet_pt);
00532   }
00533 
00534   if (jet_eta>-50.) {
00535     h_jet_eta  ->Fill(jet_eta); 
00536   }
00537 
00538   if (jet_phi>-10.) {
00539     h_jet_phi   ->Fill(jet_phi);
00540   }
00541 
00542   if (jet2_et>-10.0) {
00543     h_jet2_et  ->Fill(jet2_et);
00544   }
00545 
00546   //if (jet2_pt>0.) {
00547   //  h_jet2_pt   ->Fill(jet2_pt);
00548   //}
00549 
00550   if (jet2_eta>-50.) {
00551     h_jet2_eta  ->Fill(jet2_eta); 
00552   }
00553 
00554   if (jet2_phi>-10.) {
00555     h_jet2_phi   ->Fill(jet2_phi);
00556   }
00557  
00558 
00559 
00560   if (fill_e1 || fill_m1) {
00561     h_vertex_number->Fill(vertex_number);
00562     h_vertex_chi2->Fill(vertex_chi2);
00563     h_vertex_d0  ->Fill(vertex_d0);
00564     h_vertex_numTrks->Fill(vertex_numTrks);
00565     h_vertex_sumTrks->Fill(vertex_sumTrks);
00566   }
00567 
00568   if (fill_e1) {
00569     h_e1_et      ->Fill(electron_et);
00570     h_e1_eta     ->Fill(electron_eta);
00571     h_e1_phi     ->Fill(electron_phi);
00572   }
00573   if (fill_e2) {
00574     h_e2_et      ->Fill(electron2_et);
00575     h_e2_eta     ->Fill(electron2_eta);
00576     h_e2_phi     ->Fill(electron2_phi);
00577   }
00578   if (fill_m1) {
00579     h_m1_pt      ->Fill(muon_pt);
00580     h_m1_eta     ->Fill(muon_eta);
00581     h_m1_phi     ->Fill(muon_phi);
00582   }
00583   if (fill_m2) {
00584     h_m2_pt      ->Fill(muon2_pt);
00585     h_m2_eta     ->Fill(muon2_eta);
00586     h_m2_phi     ->Fill(muon2_phi);
00587   }
00588   if (fill_met) {
00589     h_met        ->Fill(missing_et);
00590     h_met_phi    ->Fill(met_phi);
00591   }
00593 }
00594 
00595 
00596 void EwkDQM::endJob(void) {}
00597 
00598 
00599 // This always returns only a positive deltaPhi
00600 double EwkDQM::calcDeltaPhi(double phi1, double phi2) {
00601 
00602   double deltaPhi = phi1 - phi2;
00603 
00604   if (deltaPhi < 0) deltaPhi = -deltaPhi;
00605 
00606   if (deltaPhi > 3.1415926) {
00607     deltaPhi = 2 * 3.1415926 - deltaPhi;
00608   }
00609 
00610   return deltaPhi;
00611 }
/data/refman/pasoursint/CMSSW_5_3_3/src/DQM/Physics/src/EwkDQM.cc
