Reference Manual

Go to the documentation of this file.
00001 #include "FWCore/Framework/interface/EDProducer.h"
00002 #include "FWCore/Utilities/interface/InputTag.h"
00003 #include "DataFormats/Common/interface/Handle.h"
00004 #include "FWCore/Framework/interface/Event.h"
00005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00006 #include "DataFormats/PatCandidates/interface/Muon.h"
00007 #include "FWCore/Utilities/interface/EDMException.h"
00008 #include "CommonTools/UtilAlgos/interface/TFileService.h"
00009 
00010 #include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
00011 #include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
00012 #include "DataFormats/PatCandidates/interface/Isolation.h"
00013 #include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
00014 #include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
00015 
00016 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00017 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00018 #include "DataFormats/VertexReco/interface/Vertex.h"
00019 
00020 #include <vector>
00021 
00022 using namespace edm;
00023 using namespace std;
00024 using namespace reco;
00025 using namespace isodeposit;
00026 //using namespace pat;
00027 
00028 class ZMuMuMuonUserData : public edm::EDProducer {
00029 public:
00030   ZMuMuMuonUserData( const edm::ParameterSet & );   
00031 private:
00032   void produce( edm::Event &, const edm::EventSetup & );
00033   
00034   InputTag src_,beamSpot_, primaryVertices_;
00035   double alpha_, beta_; 
00036   double ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_;
00037   string hltPath_; 
00038   template<typename T>
00039   vector<double> isolation(const T & t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta);
00040 };
00041 
00042 template<typename T>
00043 vector<double> ZMuMuMuonUserData::isolation(const T & t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta) {
00044 
00045   vector<double> iso;
00046   const pat::IsoDeposit * trkIso = t.isoDeposit(pat::TrackIso);
00047   const pat::IsoDeposit * ecalIso = t.isoDeposit(pat::EcalIso);
00048   const pat::IsoDeposit * hcalIso = t.isoDeposit(pat::HcalIso);  
00049   
00050   Direction dir = Direction(t.eta(), t.phi());
00051    
00052   pat::IsoDeposit::AbsVetos vetosTrk;
00053   vetosTrk.push_back(new ConeVeto( dir, dRVetoTrk ));
00054   vetosTrk.push_back(new ThresholdVeto( ptThreshold ));
00055   
00056   pat::IsoDeposit::AbsVetos vetosEcal;
00057   vetosEcal.push_back(new ConeVeto( dir, 0.));
00058   vetosEcal.push_back(new ThresholdVeto( etEcalThreshold ));
00059   
00060   pat::IsoDeposit::AbsVetos vetosHcal;
00061   vetosHcal.push_back(new ConeVeto( dir, 0. ));
00062   vetosHcal.push_back(new ThresholdVeto( etHcalThreshold ));
00063   
00064   double isovalueTrk = (trkIso->sumWithin(dRTrk,vetosTrk));
00065   double isovalueEcal = (ecalIso->sumWithin(dREcal,vetosEcal));
00066   double isovalueHcal = (hcalIso->sumWithin(dRHcal,vetosHcal));
00067     
00068   iso.push_back(isovalueTrk);
00069   //cout<<"isoTrk"<<isovalueTrk<<" "<<t.trackIso()<<endl;
00070   iso.push_back(isovalueEcal);
00071   //cout<<"isoEcal"<<isovalueEcal<<" "<<t.ecalIso()<<endl;
00072   iso.push_back(isovalueHcal);
00073   //cout<<"isoHcal"<<isovalueHcal<<" "<<t.hcalIso()<<endl;
00074   //double isovalueTrk  = t.trackIso();
00075   //double isovalueEcal = t.ecalIso();
00076   //double isovalueHcal = t.hcalIso();
00077 
00078   //double iso =  isovalueTrk + isovalueEcal + isovalueHcal;
00079   double combIso = alpha*( ((1+beta)/2*isovalueEcal) + ((1-beta)/2*isovalueHcal) ) + ((1-alpha)*isovalueTrk);
00080   iso.push_back(combIso);
00081   //cout<<"combIso"<<iso[3]<<endl;  
00082 
00083   double relIso = combIso /= t.pt();
00084   iso.push_back(relIso);
00085   //cout<<"relIso"<<iso[4]<<endl;
00086   return iso;
00087 }
00088 
00089 ZMuMuMuonUserData::ZMuMuMuonUserData( const ParameterSet & cfg ):
00090   src_( cfg.getParameter<InputTag>( "src" ) ),
00091   beamSpot_(cfg.getParameter<InputTag>( "beamSpot" ) ),
00092   primaryVertices_(cfg.getParameter<InputTag>( "primaryVertices" ) ),
00093   alpha_(cfg.getParameter<double>("alpha") ),
00094   beta_(cfg.getParameter<double>("beta") ), 
00095   ptThreshold_(cfg.getParameter<double >("ptThreshold") ),
00096   etEcalThreshold_(cfg.getParameter<double >("etEcalThreshold") ),
00097   etHcalThreshold_(cfg.getParameter<double >("etHcalThreshold") ),
00098   dRVetoTrk_(cfg.getParameter<double >("dRVetoTrk") ),
00099   dRTrk_(cfg.getParameter<double >("dRTrk") ),
00100   dREcal_(cfg.getParameter<double >("dREcal") ),
00101   dRHcal_(cfg.getParameter<double >("dRHcal") ),
00102   hltPath_(cfg.getParameter<std::string >("hltPath") ){
00103   produces<std::vector<pat::Muon> >();
00104 }
00105 
00106 void ZMuMuMuonUserData::produce( Event & evt, const EventSetup & ) {
00107   Handle<vector<pat::Muon>  > muons;
00108   evt.getByLabel(src_,muons);
00109 
00110   Handle<BeamSpot> beamSpotHandle;
00111   evt.getByLabel(beamSpot_, beamSpotHandle);
00112 
00113   Handle<VertexCollection> primaryVertices;  // Collection of primary Vertices
00114   evt.getByLabel(primaryVertices_, primaryVertices);
00115 
00116   auto_ptr<vector<pat::Muon> > muonColl( new vector<pat::Muon> (*muons) );
00117   for (unsigned int i = 0; i< muonColl->size();++i){
00118     pat::Muon & m = (*muonColl)[i];
00119     //pat::Muon *mu = new pat::Muon(m); 
00120     vector<double> iso = isolation(m,ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_, alpha_, beta_);
00121     m.setIsolation(pat::User1Iso, iso[0]);
00122     //cout<<"muon User1Iso "<<iso[0]<<endl;
00123     m.setIsolation(pat::User2Iso, iso[1]);
00124     //cout<<"iso2 "<<iso[1]<<endl;
00125     m.setIsolation(pat::User3Iso, iso[2]);
00126     //cout<<"iso3 "<<iso[2]<<endl;
00127     m.setIsolation(pat::User4Iso, iso[3]);
00128     //cout<<"iso4 "<<iso[3]<<endl; 
00129     m.setIsolation(pat::User5Iso, iso[4]);
00130     //cout<<"iso5 "<<iso[4]<<endl;
00131     float  zDauMuEnergyEm =  m.calEnergy().em;
00132     float  zDauMuEnergyHad = m.calEnergy().had;
00133     TrackRef muGlbRef = m.globalTrack();
00134     TrackRef muTrkRef = m.innerTrack();
00135     TrackRef muSaRef = m.outerTrack();
00136     float zDaudxyFromBS = -1;
00137     float zDaudzFromBS = -1;
00138     float zDaudxyFromPV = -1;
00139     float zDaudzFromPV = -1;
00140     int zDauNofMuChambers =   m.numberOfChambers(); 
00141     int zDauNofMuMatches = m.numberOfMatches(); 
00142     //  for the following variables looking at global/trk and sta at the same time 
00143     float zDauChi2 = -1;
00144     float zDauTrkChi2 = -1;
00145     float zDauSaChi2 = -1;
00146     float zDauNofMuonHits = -1; 
00147     float zDauSaNofMuonHits = -1; 
00148     float zDauNofStripHits = -1; 
00149     float zDauTrkNofStripHits = -1; 
00150     float zDauNofPixelHits = -1; 
00151     float zDauTrkNofPixelHits = -1; 
00152     
00153     
00154 
00155 
00156     if (muGlbRef.isNonnull() && m.isGlobalMuon() == true){
00157       zDaudxyFromBS = muGlbRef->dxy(beamSpotHandle->position());
00158       zDaudzFromBS = muGlbRef->dz(beamSpotHandle->position());
00159       zDaudxyFromPV = muGlbRef->dxy(primaryVertices->begin()->position() );
00160       zDaudzFromPV = muGlbRef->dz(primaryVertices->begin()->position() );
00161       zDauChi2 = muGlbRef->normalizedChi2();
00162       zDauTrkChi2 = muTrkRef->normalizedChi2();
00163       zDauSaChi2 = muSaRef->normalizedChi2();
00164       zDauNofMuonHits = muGlbRef->hitPattern().numberOfValidMuonHits();
00165       zDauSaNofMuonHits = muSaRef->hitPattern().numberOfValidMuonHits();
00166       zDauNofStripHits = muGlbRef->hitPattern().numberOfValidStripHits();
00167       zDauTrkNofStripHits = muTrkRef->hitPattern().numberOfValidStripHits();
00168       zDauNofPixelHits = muGlbRef->hitPattern().numberOfValidPixelHits();
00169       zDauTrkNofPixelHits = muTrkRef->hitPattern().numberOfValidPixelHits();
00170     }
00171  else if (muSaRef.isNonnull() && m.isStandAloneMuon() == true){
00172       zDaudxyFromBS = muSaRef->dxy(beamSpotHandle->position());
00173       zDaudzFromBS = muSaRef->dz(beamSpotHandle->position());
00174       zDaudxyFromPV = muSaRef->dxy(primaryVertices->begin()->position() );
00175       zDaudzFromPV = muSaRef->dz(primaryVertices->begin()->position() );
00176       zDauSaChi2 = muSaRef->normalizedChi2();
00177       zDauSaNofMuonHits = muSaRef->hitPattern().numberOfValidMuonHits(); 
00178 
00179     }  
00180      else if (muTrkRef.isNonnull() && m.isTrackerMuon() == true){
00181       zDaudxyFromBS = muTrkRef->dxy(beamSpotHandle->position());
00182       zDaudzFromBS = muTrkRef->dz(beamSpotHandle->position());
00183       zDaudxyFromPV = muTrkRef->dxy(primaryVertices->begin()->position() );
00184       zDaudzFromPV = muTrkRef->dz(primaryVertices->begin()->position() );
00185       zDauTrkChi2 = muTrkRef->normalizedChi2();
00186       zDauTrkNofStripHits = muTrkRef->hitPattern().numberOfValidStripHits();
00187       zDauTrkNofPixelHits = muTrkRef->hitPattern().numberOfValidPixelHits();
00188 
00189     }   
00190    
00191     const pat::TriggerObjectStandAloneCollection muHLTMatches =  m.triggerObjectMatchesByPath( hltPath_);
00192     float muHLTBit;
00193     int dimTrig = muHLTMatches.size();
00194         if(dimTrig !=0 ){
00195           muHLTBit = 1;
00196         } else {
00197           muHLTBit = 0; 
00198         }
00199     m.addUserFloat("zDau_dxyFromBS", zDaudxyFromBS);
00200     m.addUserFloat("zDau_dzFromBS", zDaudzFromBS);
00201     m.addUserFloat("zDau_dxyFromPV", zDaudxyFromPV);
00202     m.addUserFloat("zDau_dzFromPV", zDaudzFromPV);
00203     m.addUserFloat("zDau_HLTBit",muHLTBit);
00204     m.addUserFloat("zDau_dzFromPV", zDaudzFromPV);
00205     m.addUserFloat("zDau_Chi2", zDauChi2);
00206     m.addUserFloat("zDau_TrkChi2", zDauTrkChi2);
00207     m.addUserFloat("zDau_SaChi2", zDauSaChi2);
00208     m.addUserFloat("zDau_NofMuonHits" , zDauNofMuonHits );
00209     m.addUserFloat("zDau_SaNofMuonHits" , zDauSaNofMuonHits ); 
00210     m.addUserFloat("zDau_NofStripHits" , zDauNofStripHits ); 
00211     m.addUserFloat("zDau_TrkNofStripHits" , zDauTrkNofStripHits ); 
00212     m.addUserFloat("zDau_NofPixelHits" , zDauNofPixelHits ); 
00213     m.addUserFloat("zDau_TrkNofPixelHits" , zDauTrkNofPixelHits ); 
00214     m.addUserFloat("zDau_NofMuChambers" , zDauNofMuChambers ); 
00215     m.addUserFloat("zDau_NofMuMatches" , zDauNofMuMatches ); 
00216     m.addUserFloat("zDau_MuEnergyEm", zDauMuEnergyEm ); 
00217     m.addUserFloat("zDau_MuEnergyHad", zDauMuEnergyHad ); 
00218  }
00219 
00220   evt.put( muonColl);
00221 }
00222 
00223 #include "FWCore/Framework/interface/MakerMacros.h"
00224 
00225 DEFINE_FWK_MODULE( ZMuMuMuonUserData );
00226
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/ElectroWeakAnalysis/Skimming/plugins/ZMuMuMuonUserData.cc
