CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
TracksterP4FromEnergySumPlugin.cc
Go to the documentation of this file.
1 // Plugin for getting the four-vector from a Trackster from a simple energy sum and weighted cluster position.
2 // A simplistic 1/N(tracksters) sharing is applied for hits that belong to multiple tracksters.
3 // Alternatively takes the energy value from the pre-calculated regressed energy value in the Trackster.
4 
9 
10 namespace ticl {
12  public:
14  void setP4(const std::vector<const Trackster*>& tracksters,
15  std::vector<TICLCandidate>& ticl_cands,
16  edm::Event& event) const override;
17 
18  private:
19  std::tuple<TracksterMomentumPluginBase::LorentzVector, float> calcP4(
20  const ticl::Trackster& trackster,
21  const reco::Vertex& vertex,
22  const std::vector<reco::CaloCluster>& calo_clusters) const;
26  };
27 
29  : TracksterMomentumPluginBase(ps, std::move(ic)),
30  energy_from_regression_(ps.getParameter<bool>("energyFromRegression")),
31  vertex_token_(ic.consumes<std::vector<reco::Vertex>>(ps.getParameter<edm::InputTag>("vertices"))),
32  layer_clusters_token_(
33  ic.consumes<std::vector<reco::CaloCluster>>(ps.getParameter<edm::InputTag>("layerClusters"))) {}
34 
35  void TracksterP4FromEnergySum::setP4(const std::vector<const Trackster*>& tracksters,
36  std::vector<TICLCandidate>& ticl_cands,
37  edm::Event& event) const {
38  // Find best vertex
40  event.getByToken(vertex_token_, vertex_h);
41  auto vertex_coll = *vertex_h;
42  reco::Vertex best_vertex;
43  if (!vertex_coll.empty()) {
44  const auto& vertex = vertex_coll[0];
45  if (vertex.isValid() && !(vertex.isFake())) {
46  best_vertex = vertex;
47  }
48  }
49 
51  event.getByToken(layer_clusters_token_, layer_clusters_h);
52 
53  auto size = std::min(tracksters.size(), ticl_cands.size());
54  for (size_t i = 0; i < size; ++i) {
55  const auto* trackster = tracksters[i];
56  auto ret = calcP4(*trackster, best_vertex, *layer_clusters_h);
57 
58  auto& ticl_cand = ticl_cands[i];
59  ticl_cand.setP4(std::get<0>(ret));
60  ticl_cand.setRawEnergy(std::get<1>(ret));
61  }
62  }
63 
64  std::tuple<TracksterMomentumPluginBase::LorentzVector, float> TracksterP4FromEnergySum::calcP4(
65  const ticl::Trackster& trackster,
66  const reco::Vertex& vertex,
67  const std::vector<reco::CaloCluster>& calo_clusters) const {
68  std::array<double, 3> barycentre{{0., 0., 0.}};
69  double energy = 0.;
70  size_t counter = 0;
71 
72  for (auto idx : trackster.vertices()) {
73  auto n_vertices = trackster.vertex_multiplicity(counter++);
74  auto fraction = n_vertices ? 1.f / n_vertices : 1.f;
75  auto weight = calo_clusters[idx].energy() * fraction;
76  energy += weight;
77  barycentre[0] += calo_clusters[idx].x() * weight;
78  barycentre[1] += calo_clusters[idx].y() * weight;
79  barycentre[2] += calo_clusters[idx].z() * weight;
80  }
82  std::begin(barycentre), std::end(barycentre), std::begin(barycentre), [&energy](double val) -> double {
83  return energy >= 0. ? val / energy : val;
84  });
85 
86  math::XYZVector direction(barycentre[0] - vertex.x(), barycentre[1] - vertex.y(), barycentre[2] - vertex.z());
87  direction = direction.Unit();
88  auto raw_energy = energy;
89  energy = energy_from_regression_ ? trackster.regressed_energy() : raw_energy;
90  direction *= energy;
91 
92  math::XYZTLorentzVector cartesian(direction.X(), direction.Y(), direction.Z(), energy);
93  // Convert px, py, pz, E vector to CMS standard pt/eta/phi/m vector
95  return std::tuple(p4, raw_energy);
96  }
97 } // namespace ticl
98 
tuple ret
prodAgent to be discontinued
edm::EDGetTokenT< std::vector< reco::CaloCluster > > layer_clusters_token_
double y() const
y coordinate
Definition: Vertex.h:131
std::tuple< TracksterMomentumPluginBase::LorentzVector, float > calcP4(const ticl::Trackster &trackster, const reco::Vertex &vertex, const std::vector< reco::CaloCluster > &calo_clusters) const
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
TracksterP4FromEnergySum(const edm::ParameterSet &, edm::ConsumesCollector &&iC)
def move
Definition: eostools.py:511
double z() const
z coordinate
Definition: Vertex.h:133
T min(T a, T b)
Definition: MathUtil.h:58
double x() const
x coordinate
Definition: Vertex.h:129
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
void setP4(const std::vector< const Trackster * > &tracksters, std::vector< TICLCandidate > &ticl_cands, edm::Event &event) const override
edm::EDGetTokenT< std::vector< reco::Vertex > > vertex_token_
static std::atomic< unsigned int > counter
string end
Definition: dataset.py:937
#define DEFINE_EDM_PLUGIN(factory, type, name)
int weight
Definition: histoStyle.py:51
reco::Candidate::LorentzVector LorentzVector
tuple size
Write out results.
unsigned transform(const HcalDetId &id, unsigned transformCode)