CMS 3D CMS Logo

PrimaryVertexSorting.cc
Go to the documentation of this file.
6 #include <fastjet/internal/base.hh>
7 #include "fastjet/PseudoJet.hh"
8 #include "fastjet/JetDefinition.hh"
9 #include "fastjet/ClusterSequence.hh"
10 #include "fastjet/Selector.hh"
11 #include "fastjet/PseudoJet.hh"
13 
14 using namespace fastjet;
15 using namespace std;
16 
17 
18 float PrimaryVertexSorting::score(const reco::Vertex & pv,const std::vector<const reco::Candidate *> & cands, bool useMet ) const {
20  float sumPt2=0;
21  float sumEt=0;
22  LorentzVector met;
23  std::vector<fastjet::PseudoJet> fjInputs_;
24  fjInputs_.clear();
25  size_t countScale0 = 0;
26  for (size_t i = 0 ; i < cands.size(); i++) {
27  const reco::Candidate * c= cands[i];
28  float scale=1.;
29  if(c->bestTrack() != nullptr)
30  {
31  if(c->pt()!=0) {
32  scale=(c->pt()-c->bestTrack()->ptError())/c->pt();
33  }
34  if(edm::isNotFinite(scale)) {
35  edm::LogWarning("PrimaryVertexSorting") << "Scaling is NAN ignoring this candidate/track" << std::endl;
36  scale=0;
37  }
38  if(scale<0){
39  scale=0;
40  countScale0++;
41  }
42  }
43 
44  int absId=abs(c->pdgId());
45  if(absId==13 or absId == 11) {
46  float pt =c->pt()*scale;
47  sumPt2+=pt*pt;
48  met+=c->p4()*scale;
49  sumEt+=c->pt()*scale;
50  } else {
51  if (scale != 0){ // otherwise, what is the point to cluster zeroes
52  fjInputs_.push_back(fastjet::PseudoJet(c->px()*scale,c->py()*scale,c->pz()*scale,c->p4().E()*scale));
53  // fjInputs_.back().set_user_index(i);
54  }
55  }
56  }
57  fastjet::ClusterSequence sequence( fjInputs_, JetDefinition(antikt_algorithm, 0.4));
58  auto jets = fastjet::sorted_by_pt(sequence.inclusive_jets(0));
59  for (const auto & pj : jets) {
60  auto p4 = LorentzVector( pj.px(), pj.py(), pj.pz(), pj.e() ) ;
61  sumPt2+=(p4.pt()*p4.pt())*0.8*0.8;
62  met+=p4;
63  sumEt+=p4.pt();
64  }
65  float metAbove = met.pt() - 2*sqrt(sumEt);
66  if(metAbove > 0 and useMet) {
67  sumPt2+=metAbove*metAbove;
68  }
69  if (countScale0 == cands.size()) sumPt2 = countScale0*0.01; //leave some epsilon value to sort vertices with unknown pt
70  return sumPt2;
71 }
72 
73 
virtual double pz() const =0
z coordinate of momentum vector
virtual const Track * bestTrack() const
Definition: Candidate.h:247
constexpr bool isNotFinite(T x)
Definition: isFinite.h:9
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
virtual double py() const =0
y coordinate of momentum vector
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
virtual int pdgId() const =0
PDG identifier.
T sqrt(T t)
Definition: SSEVec.h:18
double p4[4]
Definition: TauolaWrapper.h:92
vector< PseudoJet > jets
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:814
def pv(vc)
Definition: MetAnalyzer.py:7
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
virtual double pt() const =0
transverse momentum
met
===> hadronic RAZOR
float score(const reco::Vertex &pv, const std::vector< const reco::Candidate * > &candidates, bool useMet) const
virtual double px() const =0
x coordinate of momentum vector
math::PtEtaPhiELorentzVectorF LorentzVector