CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SoftPFMuonTagInfoProducer.cc
Go to the documentation of this file.
1 // * Author: Alberto Zucchetta
2 // * Mail: a.zucchetta@cern.ch
3 // * January 16, 2015
4 
19 // Muons
25 
27 
28 // Transient Track and IP
34 #include <cmath>
35 
37  jetToken = consumes<edm::View<reco::Jet> >(conf.getParameter<edm::InputTag>("jets"));
38  muonToken = consumes<edm::View<reco::Muon> >(conf.getParameter<edm::InputTag>("muons"));
39  vertexToken = consumes<reco::VertexCollection>(conf.getParameter<edm::InputTag>("primaryVertex"));
40  pTcut = conf.getParameter<double>("muonPt");
41  SIPsigcut = conf.getParameter<double>("muonSIPsig");
42  IPsigcut = conf.getParameter<double>("filterIpsig");
43  ratio1cut = conf.getParameter<double>("filterRatio1");
44  ratio2cut = conf.getParameter<double>("filterRatio2");
45  useFilter = conf.getParameter<bool>("filterPromptMuons");
46  produces<reco::CandSoftLeptonTagInfoCollection>();
47 }
48 
50 
52  // Declare produced collection
53  std::auto_ptr<reco::CandSoftLeptonTagInfoCollection> theMuonTagInfo(new reco::CandSoftLeptonTagInfoCollection);
54 
55  // Declare and open Jet collection
56  edm::Handle<edm::View<reco::Jet> > theJetCollection;
57  iEvent.getByToken(jetToken, theJetCollection);
58 
59  // Declare Muon collection
60  edm::Handle<edm::View<reco::Muon> > theMuonCollection;
61  iEvent.getByToken(muonToken, theMuonCollection);
62 
63  // Declare and open Vertex collection
64  edm::Handle<reco::VertexCollection> theVertexCollection;
65  iEvent.getByToken(vertexToken, theVertexCollection);
66  if(!theVertexCollection.isValid() || theVertexCollection->empty()) return;
67  const reco::Vertex* vertex=&theVertexCollection->front();
68 
69  // Biult TransientTrackBuilder
71  iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder", theTrackBuilder);
72  const TransientTrackBuilder* transientTrackBuilder=theTrackBuilder.product();
73 
74  // Loop on jets
75  for(unsigned int ij=0, nj=theJetCollection->size(); ij<nj; ij++) {
76  edm::RefToBase<reco::Jet> jetRef = theJetCollection->refAt(ij);
77  // Build TagInfo object
79  tagInfo.setJetRef(jetRef);
80  // Loop on jet daughters
81  for(unsigned int id=0, nd=jetRef->numberOfDaughters(); id<nd; ++id) {
82  edm::Ptr<reco::Candidate> lepPtr = jetRef->daughterPtr(id);
83  if(std::abs(lepPtr->pdgId())!=13) continue;
84 
85  const reco::Muon* muon(NULL);
86  // Step 1: try to access the muon from reco::PFCandidate
87  const reco::PFCandidate* pfcand=dynamic_cast<const reco::PFCandidate*>(lepPtr.get());
88  if(pfcand) {
89  muon=pfcand->muonRef().get();
90  }
91  // If not PFCandidate is available, find a match looping on the muon collection
92  else {
93  for(unsigned int im=0, nm=theMuonCollection->size(); im<nm; ++im) { // --- Begin loop on muons
94  const reco::Muon* recomuon=&theMuonCollection->at(im);
95  const pat::Muon* patmuon=dynamic_cast<const pat::Muon*>(recomuon);
96  // Step 2: try a match between reco::Candidate
97  if(patmuon) {
98  if(patmuon->originalObjectRef()==lepPtr) {
99  muon=theMuonCollection->refAt(im).get();
100  break;
101  }
102  }
103  // Step 3: try a match with dR and dpT if pat::Muon casting fails
104  else {
105  if(reco::deltaR(*recomuon, *lepPtr)<0.01 && std::abs(recomuon->pt()-lepPtr->pt())/lepPtr->pt()<0.1) {
106  muon=theMuonCollection->refAt(im).get();
107  break;
108  }
109  }
110  } // --- End loop on muons
111  }
112  if(!muon || !muon::isLooseMuon(*muon) || muon->pt()<pTcut) continue;
113  reco::TrackRef trkRef( muon->innerTrack() );
114  reco::TrackBaseRef trkBaseRef( trkRef );
115  // Build Transient Track
116  reco::TransientTrack transientTrack=transientTrackBuilder->build(trkRef);
117  // Define jet and muon vectors
118  reco::Candidate::Vector jetvect(jetRef->p4().Vect()), muonvect(muon->p4().Vect());
119  // Calculate variables
120  reco::SoftLeptonProperties properties;
121  Measurement1D ip2d = IPTools::signedTransverseImpactParameter(transientTrack, GlobalVector(jetRef->px(), jetRef->py(), jetRef->pz()), *vertex).second;
122  Measurement1D ip3d = IPTools::signedImpactParameter3D(transientTrack, GlobalVector(jetRef->px(), jetRef->py(), jetRef->pz()), *vertex).second;
123  properties.sip2dsig = ip2d.significance();
124  properties.sip3dsig = ip3d.significance();
125  properties.sip2d = ip2d.value();
126  properties.sip3d = ip3d.value();
127  properties.deltaR = reco::deltaR(*jetRef, *muon);
128  properties.ptRel = ( (jetvect-muonvect).Cross(muonvect) ).R() / jetvect.R(); // | (Pj-Pu) X Pu | / | Pj |
129  float mag = muonvect.R()*jetvect.R();
130  float dot = muon->p4().Dot(jetRef->p4());
131  properties.etaRel = -log((mag - dot)/(mag + dot)) / 2.;
132  properties.ratio = muon->pt() / jetRef->pt();
133  properties.ratioRel = muon->p4().Dot(jetRef->p4()) / jetvect.Mag2();
134  properties.p0Par = boostedPPar(muon->momentum(), jetRef->momentum());
135 
136  if(std::abs(properties.sip3dsig)>SIPsigcut) continue;
137 
138  // Filter leptons from W, Z decays
139  if(useFilter && ((std::abs(properties.sip3dsig)<IPsigcut && properties.ratio>ratio1cut) || properties.ratio>ratio2cut)) continue;
140 
141  // Insert lepton properties
142  tagInfo.insert(lepPtr, properties);
143 
144  } // --- End loop on daughters
145 
146  // Fill the TagInfo collection
147  theMuonTagInfo->push_back(tagInfo);
148  } // --- End loop on jets
149 
150  // Put the TagInfo collection in the event
151  iEvent.put(theMuonTagInfo);
152 }
153 
154 
155 // compute the lepton momentum along the jet axis, in the jet rest frame
157  static const double lepton_mass = 0.00; // assume a massless (ultrarelativistic) lepton
158  static const double jet_mass = 5.279; // use B±/B0 mass as the jet rest mass [PDG 2007 updates]
159  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double> > lepton(vector.Dot(axis) / axis.r(), ROOT::Math::VectorUtil::Perp(vector, axis), 0., lepton_mass);
160  ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double> > jet( axis.r(), 0., 0., jet_mass );
161  ROOT::Math::BoostX boost( -jet.Beta() );
162  return boost(lepton).x();
163 }
164 
T getParameter(std::string const &) const
math::XYZVector Vector
point in the space
Definition: Candidate.h:43
virtual TrackRef innerTrack() const
Definition: Muon.h:48
CandidatePtr daughterPtr(size_type i) const
reference to daughter at given position
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:464
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
T const * get() const
Returns C++ pointer to the item.
Definition: Ptr.h:160
std::pair< bool, Measurement1D > signedTransverseImpactParameter(const reco::TransientTrack &track, const GlobalVector &direction, const reco::Vertex &vertex)
Definition: IPTools.cc:50
virtual Vector momentum() const
spatial momentum vector
reco::TransientTrack build(const reco::Track *p) const
std::pair< bool, Measurement1D > signedImpactParameter3D(const reco::TransientTrack &track, const GlobalVector &direction, const reco::Vertex &vertex)
Definition: IPTools.cc:71
#define NULL
Definition: scimark2.h:8
double deltaR(const T1 &t1, const T2 &t2)
Definition: deltaR.h:48
void setJetRef(const edm::Ref< T > &jetRef)
Definition: JetTagInfo.h:25
virtual double pt() const
transverse momentum
bool isLooseMuon(const reco::Muon &)
int iEvent
Definition: GenABIO.cc:230
virtual size_t numberOfDaughters() const
number of daughters
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:120
std::vector< CandSoftLeptonTagInfo > CandSoftLeptonTagInfoCollection
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::EDGetTokenT< edm::View< reco::Jet > > jetToken
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:244
const edm::Ptr< reco::Candidate > & originalObjectRef() const
reference to original object. Returns a null reference if not available
Definition: PATObject.h:498
bool isValid() const
Definition: HandleBase.h:75
virtual float boostedPPar(const math::XYZVector &, const math::XYZVector &)
tuple conf
Definition: dbtoconf.py:185
reco::MuonRef muonRef() const
Definition: PFCandidate.cc:450
void insert(const REF &lepton, const SoftLeptonProperties &properties)
double significance() const
Definition: Measurement1D.h:32
virtual double px() const
x coordinate of momentum vector
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
const T & get() const
Definition: EventSetup.h:56
T const * product() const
Definition: ESHandle.h:86
double value() const
Definition: Measurement1D.h:28
virtual double pz() const
z coordinate of momentum vector
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:39
edm::EDGetTokenT< reco::VertexCollection > vertexToken
T dot(const Basic3DVector &v) const
Scalar product, or &quot;dot&quot; product, with a vector of same type.
SoftPFMuonTagInfoProducer(const edm::ParameterSet &conf)
edm::EDGetTokenT< edm::View< reco::Muon > > muonToken
virtual void produce(edm::Event &, const edm::EventSetup &)
virtual const LorentzVector & p4() const
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
Analysis-level muon class.
Definition: Muon.h:49
virtual double py() const
y coordinate of momentum vector
Global3DVector GlobalVector
Definition: GlobalVector.h:10
tuple log
Definition: cmsBatch.py:341