CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
RECOVertex.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <sstream>
3 #include <istream>
4 #include <fstream>
5 #include <iomanip>
6 #include <cstdlib>
7 #include <cstring>
8 
9 #include "RECOVertex.h"
10 #include "HLTMessages.h"
11 
12 static const size_t kMaxVrt = 50;
13 
15  //set parameter defaults
16  _Debug = false;
17 
18  NVrtx = 0;
19  VertexCand_x = new float[kMaxVrt];
20  VertexCand_y = new float[kMaxVrt];
21  VertexCand_z = new float[kMaxVrt];
22  VertexCand_tracks = new int[kMaxVrt];
23  VertexCand_chi2 = new float[kMaxVrt];
24  VertexCand_ndof = new float[kMaxVrt];
25 }
26 
28  delete[] VertexCand_x;
29  delete[] VertexCand_y;
30  delete[] VertexCand_z;
31  delete[] VertexCand_tracks;
32  delete[] VertexCand_chi2;
33  delete[] VertexCand_ndof;
34 }
35 
37  NVrtx = 0;
38  std::memset(VertexCand_x, '\0', kMaxVrt * sizeof(float));
39  std::memset(VertexCand_y, '\0', kMaxVrt * sizeof(float));
40  std::memset(VertexCand_z, '\0', kMaxVrt * sizeof(float));
41  std::memset(VertexCand_tracks, '\0', kMaxVrt * sizeof(int));
42  std::memset(VertexCand_chi2, '\0', kMaxVrt * sizeof(float));
43  std::memset(VertexCand_ndof, '\0', kMaxVrt * sizeof(float));
44 }
45 
46 /* Setup the analysis to put the branch-variables into the tree. */
47 void RECOVertex::setup(const edm::ParameterSet& pSet, TTree* HltTree, std::string vertexType) {
48  edm::ParameterSet myHltParams = pSet.getParameter<edm::ParameterSet>("RunParameters");
49  std::vector<std::string> parameterNames = myHltParams.getParameterNames();
50 
51  for (auto& parameterName : parameterNames) {
52  if (parameterName == "Debug")
53  _Debug = myHltParams.getParameter<bool>(parameterName);
54  }
55 
56  TString br_recoNVrt = "recoNVrt";
57  br_recoNVrt.Append(vertexType);
58  HltTree->Branch(br_recoNVrt, &NVrtx, "NVrtx/I");
59 
60  TString br_recoVrtX = "recoVrtX";
61  br_recoVrtX.Append(vertexType);
62  HltTree->Branch(br_recoVrtX, VertexCand_x, "recoVrtX[NVrtx]/F");
63 
64  TString br_recoVrtY = "recoVrtY";
65  br_recoVrtY.Append(vertexType);
66  HltTree->Branch(br_recoVrtY, VertexCand_y, "recoVrtY[NVrtx]/F");
67 
68  TString br_recoVrtZ = "recoVrtZ";
69  br_recoVrtZ.Append(vertexType);
70  HltTree->Branch(br_recoVrtZ, VertexCand_z, "recoVrtZ[NVrtx]/F");
71 
72  TString br_recoVrtNtrk = "recoVrtNtrk";
73  br_recoVrtNtrk.Append(vertexType);
74  HltTree->Branch(br_recoVrtNtrk, VertexCand_tracks, "recoVrtNtrk[NVrtx]/I");
75 
76  TString br_recoVrtChi2 = "recoVrtChi2";
77  br_recoVrtChi2.Append(vertexType);
78  HltTree->Branch(br_recoVrtChi2, VertexCand_chi2, "recoVrtChi2[NVrtx]/F");
79 
80  TString br_recoVrtNdof = "recoVrtNdof";
81  br_recoVrtNdof.Append(vertexType);
82  HltTree->Branch(br_recoVrtNdof, VertexCand_ndof, "recoVrtNdof[NVrtx]/F");
83 }
84 
85 /* **Analyze the event** */
86 void RECOVertex::analyze(edm::Handle<reco::VertexCollection> recoVertexs, TTree* HltTree) {
87  // reset the tree variables
88  clear();
89 
90  if (recoVertexs.isValid()) {
91  const reco::VertexCollection* vertexs = recoVertexs.product();
92  reco::VertexCollection::const_iterator vertex_i;
93 
94  size_t size = std::min(kMaxVrt, size_t(vertexs->size()));
95  NVrtx = size;
96 
97  int nVertexCand = 0;
98  if (_Debug)
99  std::cout << "Found " << vertexs->size() << " vertices" << std::endl;
100  for (vertex_i = vertexs->begin(); vertex_i != vertexs->end(); vertex_i++) {
101  if (nVertexCand >= NVrtx)
102  break;
103  VertexCand_x[nVertexCand] = vertex_i->x();
104  VertexCand_y[nVertexCand] = vertex_i->y();
105  VertexCand_z[nVertexCand] = vertex_i->z();
106  VertexCand_tracks[nVertexCand] = vertex_i->tracksSize();
107  VertexCand_chi2[nVertexCand] = vertex_i->chi2();
108  VertexCand_ndof[nVertexCand] = vertex_i->ndof();
109  if (_Debug) {
110  std::cout << "RECOVertex -- VX, VY VZ = " << VertexCand_x[nVertexCand] << " " << VertexCand_y[nVertexCand]
111  << " " << VertexCand_z[nVertexCand] << std::endl;
112  std::cout << "RECOVertex -- Ntracks, Chi2/Dof = " << VertexCand_tracks[nVertexCand] << " "
113  << VertexCand_chi2[nVertexCand] << " / " << VertexCand_ndof[nVertexCand] << std::endl;
114  }
115  nVertexCand++;
116  }
117  }
118 }
size
Write out results.
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
T const * product() const
Definition: Handle.h:70
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
float * VertexCand_y
Definition: RECOVertex.h:33
static const size_t kMaxVrt
Definition: RECOVertex.cc:12
void setup(const edm::ParameterSet &pSet, TTree *tree, std::string vertexType)
Definition: RECOVertex.cc:47
float * VertexCand_ndof
Definition: RECOVertex.h:36
void analyze(edm::Handle< reco::VertexCollection > recoVertexs, TTree *tree)
Definition: RECOVertex.cc:86
float * VertexCand_z
Definition: RECOVertex.h:33
bool _Debug
Definition: RECOVertex.h:39
float * VertexCand_chi2
Definition: RECOVertex.h:35
int * VertexCand_tracks
Definition: RECOVertex.h:34
bool isValid() const
Definition: HandleBase.h:70
float * VertexCand_x
Definition: RECOVertex.h:33
void clear(void)
Definition: RECOVertex.cc:36
std::vector< std::string > getParameterNames() const