CMS 3D CMS Logo

VertexTableProducer.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 //
3 // Package: PhysicsTools/NanoAOD
4 // Class: VertexTableProducer
5 //
13 //
14 // Original Author: Andrea Rizzi
15 // Created: Mon, 28 Aug 2017 09:26:39 GMT
16 //
17 //
18 
19 
20 // system include files
21 #include <memory>
22 
23 // user include files
26 
29 
32 
35 
37 
44 
45 //
46 // class declaration
47 //
48 
50  public:
51  explicit VertexTableProducer(const edm::ParameterSet&);
52  ~VertexTableProducer() override;
53 
54  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
55 
56  private:
57  void beginStream(edm::StreamID) override;
58  void produce(edm::Event&, const edm::EventSetup&) override;
59  void endStream() override;
60 
61  //virtual void beginRun(edm::Run const&, edm::EventSetup const&) override;
62  //virtual void endRun(edm::Run const&, edm::EventSetup const&) override;
63  //virtual void beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
64  //virtual void endLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) override;
65 
66  // ----------member data ---------------------------
67 
77  const double dlenMin_,dlenSigMin_;
78 
79 };
80 
81 
82 
83 //
84 // constructors and destructor
85 //
87  pvs_(consumes<std::vector<reco::Vertex>>( params.getParameter<edm::InputTag>("pvSrc") )),
88  pvsScore_(consumes<edm::ValueMap<float>>( params.getParameter<edm::InputTag>("pvSrc") )),
89  svs_(consumes<edm::View<reco::VertexCompositePtrCandidate> >( params.getParameter<edm::InputTag>("svSrc") )),
90  svCut_(params.getParameter<std::string>("svCut") , true),
91  goodPvCut_(params.getParameter<std::string>("goodPvCut") , true),
92  goodPvCutString_(params.getParameter<std::string>("goodPvCut") ),
93  pvName_(params.getParameter<std::string>("pvName") ),
94  svName_(params.getParameter<std::string>("svName") ),
95  svDoc_(params.getParameter<std::string>("svDoc") ),
96  dlenMin_(params.getParameter<double>("dlenMin") ),
97  dlenSigMin_(params.getParameter<double>("dlenSigMin") )
98 
99 {
100  produces<nanoaod::FlatTable>("pv");
101  produces<nanoaod::FlatTable>("otherPVs");
102  produces<nanoaod::FlatTable>("svs");
103  produces<edm::PtrVector<reco::Candidate> >();
104 
105 }
106 
107 
109 {
110 
111  // do anything here that needs to be done at destruction time
112  // (e.g. close files, deallocate resources etc.)
113 
114 }
115 
116 
117 //
118 // member functions
119 //
120 
121 // ------------ method called to produce the data ------------
122 
123 
124 void
126 {
127  using namespace edm;
130  iEvent.getByToken(pvs_, pvsIn);
131  iEvent.getByToken(pvsScore_, pvsScoreIn);
132  auto pvTable = std::make_unique<nanoaod::FlatTable>(1,pvName_,true);
133  pvTable->addColumnValue<float>("ndof",(*pvsIn)[0].ndof(),"main primary vertex number of degree of freedom",nanoaod::FlatTable::FloatColumn,8);
134  pvTable->addColumnValue<float>("x",(*pvsIn)[0].position().x(),"main primary vertex position x coordinate",nanoaod::FlatTable::FloatColumn,10);
135  pvTable->addColumnValue<float>("y",(*pvsIn)[0].position().y(),"main primary vertex position y coordinate",nanoaod::FlatTable::FloatColumn,10);
136  pvTable->addColumnValue<float>("z",(*pvsIn)[0].position().z(),"main primary vertex position z coordinate",nanoaod::FlatTable::FloatColumn,16);
137  pvTable->addColumnValue<float>("chi2",(*pvsIn)[0].normalizedChi2(),"main primary vertex reduced chi2",nanoaod::FlatTable::FloatColumn,8);
138  int goodPVs=0;
139  for (const auto & pv : *pvsIn) if (goodPvCut_(pv)) goodPVs++;
140  pvTable->addColumnValue<int>("npvs",(*pvsIn).size(),"total number of reconstructed primary vertices",nanoaod::FlatTable::IntColumn);
141  pvTable->addColumnValue<int>("npvsGood",goodPVs,"number of good reconstructed primary vertices. selection:"+goodPvCutString_,nanoaod::FlatTable::IntColumn);
142  pvTable->addColumnValue<float>("score",(*pvsScoreIn).get(pvsIn.id(),0),"main primary vertex score, i.e. sum pt2 of clustered objects",nanoaod::FlatTable::FloatColumn,8);
143 
144  auto otherPVsTable = std::make_unique<nanoaod::FlatTable>((*pvsIn).size() >4?3:(*pvsIn).size()-1,"Other"+pvName_,false);
145  std::vector<float> pvsz;
146  for(size_t i=1;i < (*pvsIn).size() && i < 4; i++) pvsz.push_back((*pvsIn)[i-1].position().z());
147  otherPVsTable->addColumn<float>("z",pvsz,"Z position of other primary vertices, excluding the main PV",nanoaod::FlatTable::FloatColumn,8);
148 
149 
151  iEvent.getByToken(svs_, svsIn);
152  auto selCandSv = std::make_unique<PtrVector<reco::Candidate>>();
153  std::vector<float> dlen,dlenSig,pAngle,dxy,dxySig;
154  VertexDistance3D vdist;
155  VertexDistanceXY vdistXY;
156 
157  size_t i=0;
158  const auto & PV0 = pvsIn->front();
159  for (const auto & sv : *svsIn) {
160  if (svCut_(sv)) {
162  if(dl.value() > dlenMin_ and dl.significance() > dlenSigMin_){
163  dlen.push_back(dl.value());
164  dlenSig.push_back(dl.significance());
165  edm::Ptr<reco::Candidate> c = svsIn->ptrAt(i);
166  selCandSv->push_back(c);
167  double dx = (PV0.x() - sv.vx()), dy = (PV0.y() - sv.vy()), dz = (PV0.z() - sv.vz());
168  double pdotv = (dx * sv.px() + dy*sv.py() + dz*sv.pz())/sv.p()/sqrt(dx*dx + dy*dy + dz*dz);
169  pAngle.push_back(std::acos(pdotv));
170  Measurement1D d2d = vdistXY.distance(PV0, VertexState(RecoVertex::convertPos(sv.position()), RecoVertex::convertError(sv.error())));
171  dxy.push_back(d2d.value());
172  dxySig.push_back(d2d.significance());
173  }
174  }
175  i++;
176  }
177 
178 
179  auto svsTable = std::make_unique<nanoaod::FlatTable>(selCandSv->size(),svName_,false);
180  // For SV we fill from here only stuff that cannot be created with the SimpleFlatTableProducer
181  svsTable->addColumn<float>("dlen",dlen,"decay length in cm",nanoaod::FlatTable::FloatColumn,10);
182  svsTable->addColumn<float>("dlenSig",dlenSig,"decay length significance",nanoaod::FlatTable::FloatColumn, 10);
183  svsTable->addColumn<float>("dxy", dxy, "2D decay length in cm", nanoaod::FlatTable::FloatColumn, 10);
184  svsTable->addColumn<float>("dxySig", dxySig, "2D decay length significance", nanoaod::FlatTable::FloatColumn, 10);
185  svsTable->addColumn<float>("pAngle",pAngle,"pointing angle, i.e. acos(p_SV * (SV - PV)) ",nanoaod::FlatTable::FloatColumn,10);
186 
187 
188  iEvent.put(std::move(pvTable),"pv");
189  iEvent.put(std::move(otherPVsTable),"otherPVs");
190  iEvent.put(std::move(svsTable),"svs");
191  iEvent.put(std::move(selCandSv));
192 }
193 
194 // ------------ method called once each stream before processing any runs, lumis or events ------------
195 void
197 {
198 }
199 
200 // ------------ method called once each stream after processing all runs, lumis and events ------------
201 void
203 }
204 
205 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
206 void
208  //The following says we do not know what parameters are allowed so do no validation
209  // Please change this to state exactly what you do use, even if it is no parameters
211  desc.setUnknown();
212  descriptions.addDefault(desc);
213 }
214 
215 //define this as a plug-in
reco::Vertex::Point convertPos(const GlobalPoint &p)
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:137
const edm::EDGetTokenT< std::vector< reco::Vertex > > pvs_
const std::string pvName_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:579
const std::string svDoc_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
reco::Vertex::Error convertError(const GlobalError &ge)
Definition: ConvertError.h:8
const edm::EDGetTokenT< edm::ValueMap< float > > pvsScore_
Measurement1D distance(const GlobalPoint &vtx1Position, const GlobalError &vtx1PositionError, const GlobalPoint &vtx2Position, const GlobalError &vtx2PositionError) const override
const std::string goodPvCutString_
Measurement1D distance(const GlobalPoint &vtx1Position, const GlobalError &vtx1PositionError, const GlobalPoint &vtx2Position, const GlobalError &vtx2PositionError) const override
int iEvent
Definition: GenABIO.cc:230
goodPVs
Good Primary Vertex Selection.
void addDefault(ParameterSetDescription const &psetDescription)
VertexTableProducer(const edm::ParameterSet &)
T sqrt(T t)
Definition: SSEVec.h:18
const edm::EDGetTokenT< edm::View< reco::VertexCompositePtrCandidate > > svs_
def pv(vc)
Definition: MetAnalyzer.py:6
double significance() const
Definition: Measurement1D.h:32
void produce(edm::Event &, const edm::EventSetup &) override
double value() const
Definition: Measurement1D.h:28
fixed size matrix
const std::string svName_
HLT enums.
void beginStream(edm::StreamID) override
const StringCutObjectSelector< reco::Vertex > goodPvCut_
const StringCutObjectSelector< reco::Candidate > svCut_
def move(src, dest)
Definition: eostools.py:510