/ElectroWeakAnalysis/Skimming/plugins/ZMuMuTrackUserData.cc

  1#include "FWCore/Framework/interface/EDProducer.h"
  2#include "FWCore/Utilities/interface/InputTag.h"
  3#include "DataFormats/Common/interface/Handle.h"
  4#include "FWCore/Framework/interface/Event.h"
  5#include "FWCore/ParameterSet/interface/ParameterSet.h"
  6#include "DataFormats/PatCandidates/interface/GenericParticle.h"
  7#include "FWCore/Utilities/interface/EDMException.h"
  8
  9#include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
 10#include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
 11#include "DataFormats/PatCandidates/interface/Isolation.h"
 12#include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
 13#include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
 14
 15#include "DataFormats/BeamSpot/interface/BeamSpot.h"
 16#include "DataFormats/VertexReco/interface/VertexFwd.h"
 17#include "DataFormats/VertexReco/interface/Vertex.h"
 18
 19#include <vector>
 20
 21using namespace edm;
 22using namespace std;
 23using namespace reco;
 24using namespace isodeposit;
 25//using namespace pat;
 26
 27class ZMuMuTrackUserData : public edm::EDProducer {
 28public:
 29  ZMuMuTrackUserData( const edm::ParameterSet & );   
 30private:
 31  void produce( edm::Event &, const edm::EventSetup & ) override;
 32  
 33  InputTag src_,beamSpot_, primaryVertices_;
 34  double ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_;
 35  double alpha_, beta_; 
 36  template<typename T>
 37  vector<double> isolation(const T & t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta);
 38
 39};
 40
 41template<typename T>
 42vector<double> ZMuMuTrackUserData::isolation(const T & t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta) {
 43
 44  vector<double> iso;
 45  const pat::IsoDeposit * trkIso = t.isoDeposit(pat::TrackIso);
 46  const pat::IsoDeposit * ecalIso = t.isoDeposit(pat::EcalIso);
 47  const pat::IsoDeposit * hcalIso = t.isoDeposit(pat::HcalIso);  
 48  
 49  Direction dir = Direction(t.eta(), t.phi());
 50
 51   
 52  pat::IsoDeposit::AbsVetos vetosTrk;
 53  vetosTrk.push_back(new ConeVeto( dir, dRVetoTrk ));
 54  vetosTrk.push_back(new ThresholdVeto( ptThreshold ));
 55  
 56  pat::IsoDeposit::AbsVetos vetosEcal;
 57  vetosEcal.push_back(new ConeVeto( dir, 0.));
 58  vetosEcal.push_back(new ThresholdVeto( etEcalThreshold ));
 59  
 60  pat::IsoDeposit::AbsVetos vetosHcal;
 61  vetosHcal.push_back(new ConeVeto( dir, 0. ));
 62  vetosHcal.push_back(new ThresholdVeto( etHcalThreshold ));
 63  
 64  double isovalueTrk = (trkIso->sumWithin(dRTrk,vetosTrk));
 65  double isovalueEcal = (ecalIso->sumWithin(dREcal,vetosEcal));
 66  double isovalueHcal = (hcalIso->sumWithin(dRHcal,vetosHcal));
 67
 68  iso.push_back(isovalueTrk);
 69  iso.push_back(isovalueEcal);
 70  iso.push_back(isovalueHcal);
 71
 72  //double iso =  isovalueTrk + isovalueEcal + isovalueHcal;
 73  double combIso = alpha*( ((1+beta)/2*isovalueEcal) + ((1-beta)/2*isovalueHcal) ) + ((1-alpha)*isovalueTrk);
 74 
 75  iso.push_back(combIso);
 76  double relIso = combIso /= t.pt();
 77  iso.push_back(relIso);
 78  return iso;
 79}
 80
 81ZMuMuTrackUserData::ZMuMuTrackUserData( const ParameterSet & cfg ):
 82  src_( cfg.getParameter<InputTag>( "src" ) ),
 83  beamSpot_(cfg.getParameter<InputTag>( "beamSpot" ) ),
 84  primaryVertices_(cfg.getParameter<InputTag>( "primaryVertices" ) ),
 85  ptThreshold_(cfg.getParameter<double >("ptThreshold") ),
 86  etEcalThreshold_(cfg.getParameter<double >("etEcalThreshold") ),
 87  etHcalThreshold_(cfg.getParameter<double >("etHcalThreshold") ),
 88  dRVetoTrk_(cfg.getParameter<double >("dRVetoTrk") ),
 89  dRTrk_(cfg.getParameter<double >("dRTrk") ),
 90  dREcal_(cfg.getParameter<double >("dREcal") ),
 91  dRHcal_(cfg.getParameter<double >("dRHcal") ),
 92  alpha_(cfg.getParameter<double>("alpha") ),
 93  beta_(cfg.getParameter<double>("beta") ){
 94  produces<std::vector<pat::GenericParticle> >();
 95}
 96
 97void ZMuMuTrackUserData::produce( Event & evt, const EventSetup & ) {
 98  Handle<vector<pat::GenericParticle>  > tracks;
 99  evt.getByLabel(src_,tracks);
100
101  Handle<BeamSpot> beamSpotHandle;
102  evt.getByLabel(beamSpot_, beamSpotHandle);
103  
104
105  Handle<VertexCollection> primaryVertices;  // Collection of primary Vertices
106  evt.getByLabel(primaryVertices_, primaryVertices);
107  
108  auto_ptr<vector<pat::GenericParticle> > tkColl( new vector<pat::GenericParticle> (*tracks) );
109  for (unsigned int i = 0; i< tkColl->size();++i){
110    pat::GenericParticle & tk = (*tkColl)[i];
111    vector<double> iso = isolation(tk,ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_, alpha_, beta_);
112    tk.setIsolation(pat::User1Iso, iso[0]);
113    //    cout << "track User1Iso " << iso[0] << endl; 
114    tk.setIsolation(pat::User2Iso, iso[1]);
115    //cout << "track User2Iso " << iso[1] << endl; 
116    tk.setIsolation(pat::User3Iso, iso[2]);
117    //cout << "track User3Iso " << iso[2] << endl; 
118    tk.setIsolation(pat::User4Iso, iso[3]);
119    //cout << "track User4Iso " << iso[3] << endl; 
120    tk.setIsolation(pat::User5Iso, iso[4]);
121    //cout << "track User5Iso " << iso[4] << endl; 
122  
123    
124
125    float zDaudxyFromBS = -1 ;
126    float zDaudzFromBS = -1;
127    float zDaudxyFromPV = -1;
128    float zDaudzFromPV = -1;	
129    float zDauNofMuChambers = -1;
130    float zDauNofMuMatches = -1;
131    float zDauChi2 = -1;
132    float zDauTrkChi2 = -1;
133    float zDauSaChi2 = -1;
134    float zDauNofMuonHits =- 1;
135    float zDauNofStripHits = -1;
136    float zDauNofPixelHits = -1;
137    float zDauMuEnergyEm = -1;
138    float zDauMuEnergyHad = -1;
139
140    TrackRef muTrkRef = tk.track();
141    if (muTrkRef.isNonnull()){ 
142      zDaudxyFromBS = muTrkRef->dxy(beamSpotHandle->position());
143      zDaudzFromBS = muTrkRef->dz(beamSpotHandle->position());
144      zDaudxyFromPV = muTrkRef->dxy(primaryVertices->begin()->position() );
145      zDaudzFromPV = muTrkRef->dz(primaryVertices->begin()->position() );	
146      zDauChi2 = muTrkRef->normalizedChi2();
147      zDauTrkChi2 = muTrkRef->normalizedChi2();
148      zDauNofStripHits = muTrkRef->hitPattern().numberOfValidStripHits();
149      zDauNofPixelHits = muTrkRef->hitPattern().numberOfValidPixelHits();
150    }
151    tk.addUserFloat("zDau_dxyFromBS", zDaudxyFromBS);
152    tk.addUserFloat("zDau_dzFromBS", zDaudzFromBS);
153    tk.addUserFloat("zDau_dxyFromPV", zDaudxyFromPV);
154    tk.addUserFloat("zDau_dzFromPV", zDaudzFromPV);
155    tk.addUserFloat("zDau_NofMuonHits" , zDauNofMuonHits ); 
156    tk.addUserFloat("zDau_TrkNofStripHits" , zDauNofStripHits ); 
157    tk.addUserFloat("zDau_TrkNofPixelHits" , zDauNofPixelHits ); 
158    tk.addUserFloat("zDau_NofMuChambers", zDauNofMuChambers);
159    tk.addUserFloat("zDau_NofMuMatches", zDauNofMuMatches);
160    tk.addUserFloat("zDau_Chi2", zDauChi2);
161    tk.addUserFloat("zDau_TrkChi2", zDauTrkChi2);
162    tk.addUserFloat("zDau_SaChi2", zDauSaChi2);
163    tk.addUserFloat("zDau_MuEnergyEm", zDauMuEnergyEm);
164    tk.addUserFloat("zDau_MuEnergyHad", zDauMuEnergyHad);
165
166
167  }
168
169  evt.put( tkColl);
170}
171
172#include "FWCore/Framework/interface/MakerMacros.h"
173
174DEFINE_FWK_MODULE( ZMuMuTrackUserData );
175