File indexing completed on 2024-10-29 06:08:14
0001 #ifndef ParticleFlowReco_PFRecHit_h
0002 #define ParticleFlowReco_PFRecHit_h
0003
0004
0005 #include <vector>
0006 #include <map>
0007 #include <memory>
0008 #include <iostream>
0009
0010 #include "DataFormats/Math/interface/Point3D.h"
0011 #include "DataFormats/Math/interface/Vector3D.h"
0012 #include "Math/GenVector/PositionVector3D.h"
0013 #include "DataFormats/ParticleFlowReco/interface/PFLayer.h"
0014 #include "DataFormats/ParticleFlowReco/interface/PFRecHitFwd.h"
0015
0016 #include "DataFormats/CaloRecHit/interface/CaloRecHit.h"
0017 #include "DataFormats/Common/interface/RefToBase.h"
0018
0019 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0020 #include "Geometry/CaloGeometry/interface/CaloCellGeometryMayOwnPtr.h"
0021
0022 namespace reco {
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032 class PFRecHit {
0033 public:
0034 using PositionType = GlobalPoint::BasicVectorType;
0035 using REPPoint = RhoEtaPhi;
0036 using RepCorners = CaloCellGeometry::RepCorners;
0037 using REPPointVector = RepCorners;
0038 using CornersVec = CaloCellGeometry::CornersVec;
0039
0040 struct Neighbours {
0041 using Pointer = unsigned int const*;
0042 Neighbours() {}
0043 Neighbours(Pointer ib, unsigned int n) : b(ib), e(ib + n) {}
0044 Pointer b, e;
0045 Pointer begin() const { return b; }
0046 Pointer end() const { return e; }
0047 unsigned int size() const { return e - b; }
0048 };
0049
0050 enum { NONE = 0 };
0051
0052 PFRecHit() {}
0053
0054 PFRecHit(
0055 CaloCellGeometryMayOwnPtr caloCell, unsigned int detId, PFLayer::Layer layer, float energy, uint32_t flags = 0)
0056 : caloCell_(std::move(caloCell)), detId_(detId), layer_(layer), energy_(energy), flags_(flags) {}
0057
0058
0059 PFRecHit(const PFRecHit& other) = default;
0060 PFRecHit(PFRecHit&& other) = default;
0061 PFRecHit& operator=(const PFRecHit& other) = default;
0062 PFRecHit& operator=(PFRecHit&& other) = default;
0063
0064
0065 ~PFRecHit() = default;
0066
0067 void setEnergy(float energy) { energy_ = energy; }
0068
0069 void addNeighbour(short x, short y, short z, unsigned int);
0070 unsigned int getNeighbour(short x, short y, short z) const;
0071 void setTime(double time) { time_ = time; }
0072 void setDepth(int depth) { depth_ = depth; }
0073 void clearNeighbours() {
0074 neighbours_.clear();
0075 neighbourInfos_.clear();
0076 neighbours4_ = neighbours8_ = 0;
0077 }
0078
0079 Neighbours neighbours4() const { return buildNeighbours(neighbours4_); }
0080 Neighbours neighbours8() const { return buildNeighbours(neighbours8_); }
0081
0082 Neighbours neighbours() const { return buildNeighbours(neighbours_.size()); }
0083
0084 const std::vector<unsigned short>& neighbourInfos() { return neighbourInfos_; }
0085
0086
0087 CaloCellGeometry const& caloCell() const { return *(caloCell_.get()); }
0088 bool hasCaloCell() const { return (caloCell_ != nullptr); }
0089
0090
0091 unsigned detId() const { return detId_; }
0092
0093
0094 PFLayer::Layer layer() const { return layer_; }
0095
0096
0097 float energy() const { return energy_; }
0098
0099
0100 float time() const { return time_; }
0101
0102
0103 int depth() const { return depth_; }
0104
0105
0106 double pt2() const { return energy_ * energy_ * (position().perp2() / position().mag2()); }
0107
0108
0109 uint32_t flags() const { return flags_; }
0110
0111
0112 void setFlags(uint32_t flags) { flags_ = flags; }
0113
0114
0115 PositionType const& position() const { return caloCell().getPosition().basicVector(); }
0116
0117 RhoEtaPhi const& positionREP() const { return caloCell().repPos(); }
0118
0119
0120 CornersVec const& getCornersXYZ() const { return caloCell().getCorners(); }
0121
0122 RepCorners const& getCornersREP() const { return caloCell().getCornersREP(); }
0123
0124
0125 bool operator>=(const PFRecHit& rhs) const { return (energy_ >= rhs.energy_); }
0126
0127
0128 bool operator>(const PFRecHit& rhs) const { return (energy_ > rhs.energy_); }
0129
0130
0131 bool operator<=(const PFRecHit& rhs) const { return (energy_ <= rhs.energy_); }
0132
0133
0134 bool operator<(const PFRecHit& rhs) const { return (energy_ < rhs.energy_); }
0135
0136 private:
0137 Neighbours buildNeighbours(unsigned int n) const { return Neighbours(neighbours_.data(), n); }
0138
0139
0140 CaloCellGeometryMayOwnPtr caloCell_;
0141
0142
0143 unsigned int detId_ = 0;
0144
0145
0146 PFLayer::Layer layer_ = PFLayer::NONE;
0147
0148
0149 float energy_ = 0;
0150
0151
0152 float time_ = -1;
0153
0154
0155 int depth_ = 0;
0156
0157
0158 std::vector<unsigned int> neighbours_;
0159 std::vector<unsigned short> neighbourInfos_;
0160
0161
0162 unsigned int neighbours4_ = 0;
0163 unsigned int neighbours8_ = 0;
0164
0165
0166 uint32_t flags_ = 0;
0167 };
0168
0169 }
0170 std::ostream& operator<<(std::ostream& out, const reco::PFRecHit& hit);
0171
0172 #endif