CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch13/src/HLTriggerOffline/SUSYBSM/src/McSelector.cc

Go to the documentation of this file.
00001 /*  \class McSelector
00002  *
00003  *  Class to apply analysis cuts in the TriggerValidation Code
00004  *
00005  *  Author: Massimiliano Chiorboli      Date: August 2007
00006  //         Maurizio Pierini
00007  //         Maria Spiropulu
00008  //    Philip Hebda, July 2009
00009 *
00010  */
00011 
00012 #include "HLTriggerOffline/SUSYBSM/interface/McSelector.h"
00013 
00014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00015 
00016 using namespace edm;
00017 using namespace reco;
00018 using namespace std;
00019 
00020 McSelector::McSelector(edm::ParameterSet userCut_params)
00021 {
00022   //******************** PLEASE PAY ATTENTION: number of electron and muons is strictly equal, for jets, taus and photons the requirement is >= ********************
00023   name     = userCut_params.getParameter<string>("name");
00024   m_genSrc       = userCut_params.getParameter<string>("mcparticles");
00025   m_genJetSrc    = userCut_params.getParameter<string>("genJets");
00026   m_genMetSrc    = userCut_params.getParameter<string>("genMet");
00027   mc_ptElecMin     = userCut_params.getParameter<double>("mc_ptElecMin"     );
00028   mc_ptMuonMin     = userCut_params.getParameter<double>("mc_ptMuonMin"     );
00029   mc_ptTauMin      = userCut_params.getParameter<double>("mc_ptTauMin"      );
00030   mc_ptPhotMin     = userCut_params.getParameter<double>("mc_ptPhotMin"     );
00031   mc_ptJetMin      = userCut_params.getParameter<double>("mc_ptJetMin"      );   
00032   mc_ptJetForHtMin = userCut_params.getParameter<double>("mc_ptJetForHtMin" );
00033   mc_metMin        = userCut_params.getParameter<double>("mc_metMin"        );
00034   mc_htMin         = userCut_params.getParameter<double>("mc_htMin"         );
00035   mc_nElec     = userCut_params.getParameter<int>("mc_nElec");
00036   mc_nElecRule = userCut_params.getParameter<string>("mc_nElecRule");
00037   mc_nMuon     = userCut_params.getParameter<int>("mc_nMuon");
00038   mc_nMuonRule = userCut_params.getParameter<string>("mc_nMuonRule");
00039   mc_nTau      = userCut_params.getParameter<int>("mc_nTau");
00040   mc_nPhot     = userCut_params.getParameter<int>("mc_nPhot");
00041   mc_nJet      = userCut_params.getParameter<int>("mc_nJet" );
00042 
00043 
00044   edm::LogInfo("HLTriggerOfflineSUSYBSM") << "UserAnalysis parameters, MC for " << name << " selection:" ;
00045   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptElecMin      = " << mc_ptElecMin    ;
00046   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptMuonMin      = " << mc_ptMuonMin    ;
00047   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptTauMin       = " << mc_ptTauMin     ;
00048   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptPhotMin      = " << mc_ptPhotMin    ;
00049   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptJetMin       = " << mc_ptJetMin     ;
00050   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_ptJetForHtMin  = " << mc_ptJetForHtMin   ;
00051   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_metMin         = " << mc_metMin       ;
00052   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_htMin          = " << mc_htMin        ;
00053 
00054   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nElec          = " << mc_nElec   ;
00055   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nElecRule = " << mc_nElecRule ;
00056   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nMuon          = " << mc_nMuon   ;
00057   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nMuonRule = " << mc_nMuonRule ;
00058   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nTau           = " << mc_nTau    ;
00059   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nPhot          = " << mc_nPhot   ;
00060   edm::LogInfo("HLTriggerOfflineSUSYBSM") << " mc_nJet           = " << mc_nJet    ;
00061 
00062 
00063 }
00064 
00065 string McSelector::GetName(){return name;}
00066 
00067 bool McSelector::isSelected(const edm::Event& iEvent)
00068 {
00069 
00070 
00071   this->handleObjects(iEvent);
00072 
00073 
00074   bool TotalCutPassed = false;
00075 
00076 
00077   bool ElectronCutPassed = false;
00078   bool MuonCutPassed = false;
00079   bool TauCutPassed = false;
00080   bool PhotonCutPassed = false;
00081   bool JetCutPassed = false;  
00082   bool MetCutPassed = false;  
00083   bool HtCutPassed = false;  
00084 
00085   int nMcElectrons = 0;
00086   int nMcMuons     = 0;
00087   int nMcTaus     = 0;
00088   int nMcPhotons   = 0;
00089   int nMcJets      = 0;
00090   
00091   //  cout <<"----------------------------------------------------------------------------" << endl;
00092   
00093   for(unsigned int i=0; i<theGenParticleCollection->size(); i++) {
00094     const GenParticle* genParticle = (&(*theGenParticleCollection)[i]);    
00095     if(genParticle->status() == 2) {
00096       //taus
00097       if(fabs(genParticle->pdgId()) == 15) {
00098         //      cout << "Tau Mother = " << genParticle->mother()->pdgId() << endl;
00099         //      if(fabs(genParticle->mother()->pdgId()) == 15)  cout << "Tau GrandMother = " << genParticle->mother()->mother()->pdgId() << endl;
00100         if(genParticle->pt() > mc_ptTauMin && fabs(genParticle->eta())<2.5) {
00101           nMcTaus++;
00102         }
00103       }
00104     }
00105     if(genParticle->status() == 1) {
00106       //      cout << "genParticle->status() = " << genParticle->status() << "    genParticle->pdgId() = " << genParticle->pdgId() << "     genParticle->pt() = " << genParticle->pt() << endl;
00107       //electrons
00108       if(fabs(genParticle->pdgId()) == 11 && genParticle->numberOfMothers()) {
00109         //              cout << "Mc Electron, pt = " << genParticle->pt() << endl;
00110         //              cout << "Electron Mother = " << genParticle->mother()->pdgId() << endl;
00111 //              if(fabs(genParticle->mother()->pdgId()) == 11 || fabs(genParticle->mother()->pdgId()) == 15) 
00112 //                cout << "Electron GrandMother = " << genParticle->mother()->mother()->pdgId() << endl;
00113 //              if(fabs(genParticle->mother()->mother()->pdgId()) == 11 || fabs(genParticle->mother()->mother()->pdgId()) == 15)
00114 //                cout << "Electron GrandGrandMother = " << genParticle->mother()->mother()->mother()->pdgId() << endl;
00115         int motherCode = genParticle->mother()->pdgId();
00116         if(fabs(motherCode) == fabs(genParticle->pdgId()) || fabs(motherCode) == 15) {
00117           motherCode = genParticle->mother()->mother()->pdgId();
00118           if(fabs(motherCode) == fabs(genParticle->pdgId()) || fabs(motherCode) == 15) motherCode = genParticle->mother()->mother()->mother()->pdgId();
00119         }
00120         if(fabs(motherCode) == 23 || fabs(motherCode) == 24 || fabs(motherCode) > 999999) {
00121           if(genParticle->pt() > mc_ptElecMin && fabs(genParticle->eta())<2.5) {
00122             nMcElectrons++;
00123             //            cout << "Mc Electron, pt = " << genParticle->pt() << endl;
00124           }
00125         }
00126       }
00127       //muons
00128       if(fabs(genParticle->pdgId()) == 13 && genParticle->numberOfMothers()) {
00129         //      cout << "Mc Muon, pt = " << genParticle->pt() << endl;
00130         //              cout << "Muon Mother = " << genParticle->mother()->pdgId() << endl;
00131 //              if(fabs(genParticle->mother()->pdgId()) == 13 || fabs(genParticle->mother()->pdgId()) == 15)
00132 //                cout << "Muon GrandMother = " << genParticle->mother()->mother()->pdgId() << endl;
00133 //              if(fabs(genParticle->mother()->mother()->pdgId()) == 13 || fabs(genParticle->mother()->mother()->pdgId()) == 15)
00134 //                cout << "Muon GrandGrandMother = " << genParticle->mother()->mother()->mother()->pdgId() << endl;
00135         int motherCode = genParticle->mother()->pdgId();
00136         if(fabs(motherCode) == fabs(genParticle->pdgId()) || fabs(motherCode) == 15) {
00137           motherCode = genParticle->mother()->mother()->pdgId();
00138           if(fabs(motherCode) == fabs(genParticle->pdgId()) || fabs(motherCode) == 15) motherCode = genParticle->mother()->mother()->mother()->pdgId();
00139         }
00140         if(fabs(motherCode) == 23 || fabs(motherCode) == 24 || fabs(motherCode) > 999999) {
00141           if(genParticle->pt() > mc_ptMuonMin && fabs(genParticle->eta())<2.5) {
00142             nMcMuons++;
00143           //              cout << "Mc Muon, pt = " << genParticle->pt() << endl;
00144           }
00145         }
00146       }
00147       //photons
00148       if(fabs(genParticle->pdgId()) == 22 && fabs(genParticle->eta())<2.5) {
00149         //      cout << "Mc Photon, pt = " << genParticle->pt() << endl;
00150         if(genParticle->pt() > mc_ptPhotMin) {
00151           //      cout << "Photon Mother = " << genParticle->mother()->pdgId() << endl;
00152           nMcPhotons++;
00153           //      cout << "Mc Photon, pt = " << genParticle->pt() << endl;
00154         }
00155       }
00156     }
00157   }
00158 
00159   ht = 0;
00160   for(unsigned int i=0; i<theGenJetCollection->size();i++) {
00161     GenJet genjet = (*theGenJetCollection)[i];
00162     //    cout << "Mc Jet, pt = " << genjet.pt() << endl;
00163     if(genjet.pt() > mc_ptJetMin && fabs(genjet.eta())<3) {
00164       nMcJets++;
00165       //           cout << "Mc Jet, pt = " << genjet.pt() << endl;
00166     }
00167     if(genjet.pt() > mc_ptJetForHtMin) ht =+ genjet.pt();
00168   }
00169   if(ht>mc_htMin) HtCutPassed = true;
00170   
00171   
00172   const GenMET theGenMET = theGenMETCollection->front();
00173   //  cout << "GenMET = " << theGenMET.pt() << endl;
00174   if(theGenMET.pt() > mc_metMin) MetCutPassed = true;
00175   
00176   //  cout <<"----------------------------------------------------------------------------" << endl;
00177   
00178 
00179 //   cout << "nMcElectrons = " << nMcElectrons << endl;
00180 //   cout << "nMcMuons     = " << nMcMuons     << endl;
00181 //   cout << "nMcTaus      = " << nMcTaus      << endl;
00182 //   cout << "nMcPhotons   = " << nMcPhotons   << endl;
00183 //   cout << "nMcJets      = " << nMcJets      << endl;
00184   
00185   if(mc_nElecRule == "strictEqual") {
00186     if(nMcElectrons == mc_nElec) ElectronCutPassed = true;
00187   } 
00188   else if (mc_nElecRule == "greaterEqual") {
00189     if(nMcElectrons >= mc_nElec) ElectronCutPassed = true;
00190   } 
00191   else cout << "WARNING: not a correct definition of cuts at gen level for electrons! " << endl;
00192   
00193   
00194   if(mc_nMuonRule == "strictEqual") {
00195     if(nMcMuons     == mc_nMuon) MuonCutPassed     = true;
00196   }
00197   else if (mc_nMuonRule == "greaterEqual") {
00198     if(nMcMuons >= mc_nMuon) MuonCutPassed = true;
00199   } 
00200   else cout << "WARNING: not a correct definition of cuts at gen level for muons! " << endl;
00201   if(nMcTaus      >= mc_nTau)  TauCutPassed      = true;
00202   if(nMcPhotons   >= mc_nPhot) PhotonCutPassed   = true;
00203   if(nMcJets      >= mc_nJet ) JetCutPassed      = true;
00204   
00205   
00206 //   cout << "ElectronCutPassed = " << (int)ElectronCutPassed << endl;
00207 //   cout << "MuonCutPassed     = " << (int)MuonCutPassed << endl;
00208 //   cout << "PhotonCutPassed   = " << (int)PhotonCutPassed << endl;
00209 //   cout << "JetCutPassed      = " << (int)JetCutPassed << endl;
00210 //   cout << "MetCutPassed      = " << (int)MetCutPassed << endl;
00211   
00212   
00213 //   if(
00214 //      (ElectronCutPassed || MuonCutPassed) &&
00215 //      PhotonCutPassed                      &&
00216 //      JetCutPassed                         &&
00217 //      MetCutPassed      )   TotalCutPassed = true;
00218   
00219   
00220   if(
00221      ElectronCutPassed   && 
00222      MuonCutPassed       &&
00223      TauCutPassed        &&
00224      PhotonCutPassed     &&
00225      JetCutPassed        &&
00226      MetCutPassed        &&
00227      HtCutPassed          )   TotalCutPassed = true;
00228 
00229 
00230   // Apply a veto: removed because we require the exact number of leptons
00231   // and not >=
00232   // the veto is hence equivalent to request 0 leptons
00233 
00234 //   if(TotalCutPassed) {
00235 //     if(mc_ElecVeto) {
00236 //       if(nMcElectrons>0) TotalCutPassed = false;
00237 //     }
00238 //     if(mc_MuonVeto) {
00239 //       if(nMcMuons>0) TotalCutPassed = false;
00240 //     }
00241 //   }
00242 
00243   
00244 //  cout << "TotalCutPassed = " << TotalCutPassed << endl;
00245 
00246   return TotalCutPassed;
00247 }
00248  
00249 void McSelector::handleObjects(const edm::Event& iEvent)
00250 {
00251 
00252   //Get the GenParticleCandidates
00253   Handle< reco::GenParticleCollection > theCandidateCollectionHandle;
00254   iEvent.getByLabel(m_genSrc, theCandidateCollectionHandle);
00255   theGenParticleCollection = theCandidateCollectionHandle.product();
00256 
00257   
00258   //Get the GenJets
00259   Handle< GenJetCollection > theGenJetCollectionHandle ;
00260   iEvent.getByLabel( m_genJetSrc, theGenJetCollectionHandle);  
00261   theGenJetCollection = theGenJetCollectionHandle.product();
00262 
00263 
00264   //Get the GenMET
00265   Handle< GenMETCollection > theGenMETCollectionHandle;
00266   iEvent.getByLabel( m_genMetSrc, theGenMETCollectionHandle);
00267   theGenMETCollection = theGenMETCollectionHandle.product();
00268 
00269 
00270 }