CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
Cluster3DPCACalculator Class Reference
Inheritance diagram for Cluster3DPCACalculator:
PFCPositionCalculatorBase

Public Member Functions

void calculateAndSetPosition (reco::PFCluster &, const HcalPFCuts *) override
 
void calculateAndSetPositions (reco::PFClusterCollection &, const HcalPFCuts *) override
 
 Cluster3DPCACalculator (const edm::ParameterSet &conf, edm::ConsumesCollector &cc)
 
 Cluster3DPCACalculator (const Cluster3DPCACalculator &)=delete
 
Cluster3DPCACalculatoroperator= (const Cluster3DPCACalculator &)=delete
 
- Public Member Functions inherited from PFCPositionCalculatorBase
const std::string & name () const
 
PosCalcoperator= (const PosCalc &)=delete
 
 PFCPositionCalculatorBase (const edm::ParameterSet &conf, edm::ConsumesCollector &cc)
 
 PFCPositionCalculatorBase (const PosCalc &)=delete
 
virtual void update (const edm::EventSetup &)
 
virtual ~PFCPositionCalculatorBase ()=default
 

Private Member Functions

void calculateAndSetPositionActual (reco::PFCluster &)
 
void showerParameters (const reco::PFCluster &, math::XYZPoint &, math::XYZVector &)
 

Private Attributes

std::unique_ptr< TPrincipal > pca_
 
const bool updateTiming_
 

Additional Inherited Members

- Protected Attributes inherited from PFCPositionCalculatorBase
const float _minFractionInCalc
 

Detailed Description

Definition at line 21 of file Cluster3DPCACalculator.cc.

Constructor & Destructor Documentation

◆ Cluster3DPCACalculator() [1/2]

Cluster3DPCACalculator::Cluster3DPCACalculator ( const edm::ParameterSet conf,
edm::ConsumesCollector cc 
)
inline

Definition at line 23 of file Cluster3DPCACalculator.cc.

25  updateTiming_(conf.getParameter<bool>("updateTiming")),
26  pca_(new TPrincipal(3, "D")) {}
PFCPositionCalculatorBase(const edm::ParameterSet &conf, edm::ConsumesCollector &cc)
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
uint32_t cc[maxCellsPerHit]
Definition: gpuFishbone.h:49
std::unique_ptr< TPrincipal > pca_

◆ Cluster3DPCACalculator() [2/2]

Cluster3DPCACalculator::Cluster3DPCACalculator ( const Cluster3DPCACalculator )
delete

Member Function Documentation

◆ calculateAndSetPosition()

void Cluster3DPCACalculator::calculateAndSetPosition ( reco::PFCluster cluster,
const HcalPFCuts cuts 
)
overridevirtual

Implements PFCPositionCalculatorBase.

Definition at line 44 of file Cluster3DPCACalculator.cc.

References calculateAndSetPositionActual(), and pca_.

44  {
45  pca_ = std::make_unique<TPrincipal>(3, "D");
47 }
void calculateAndSetPositionActual(reco::PFCluster &)
std::unique_ptr< TPrincipal > pca_

◆ calculateAndSetPositionActual()

void Cluster3DPCACalculator::calculateAndSetPositionActual ( reco::PFCluster cluster)
private

Definition at line 56 of file Cluster3DPCACalculator.cc.

References reco::PFCluster::calculatePositionREP(), Exception, mps_fire::i, createfilelist::int, edm::isNotFinite(), compareTotals::means, SiStripPI::min, PFLayer::NONE, pca_, reco::PFCluster::recHitFractions(), reco::CaloCluster::seed(), reco::CaloCluster::setEnergy(), reco::PFCluster::setLayer(), reco::CaloCluster::setPosition(), reco::PFCluster::setTime(), and updateTiming_.

Referenced by calculateAndSetPosition(), and calculateAndSetPositions().

56  {
57  if (!cluster.seed()) {
58  throw cms::Exception("ClusterWithNoSeed") << " Found a cluster with no seed: " << cluster;
59  }
60  double cl_energy = 0;
61  double max_e = 0.0;
62  double avg_time = 0.0;
63  double time_norm = 0.0;
64  PFLayer::Layer max_e_layer = PFLayer::NONE;
65  reco::PFRecHitRef refseed;
66  double pcavars[3];
67 
68  for (const reco::PFRecHitFraction& rhf : cluster.recHitFractions()) {
69  const reco::PFRecHitRef& refhit = rhf.recHitRef();
70  double rh_energy = refhit->energy();
71  double rh_time = refhit->time();
72  cl_energy += rh_energy * rhf.fraction();
73  if (rh_time > 0.0) { // time == -1 means no measurement
74  // all times are offset by one nanosecond in digitizer
75  // remove that here so all times of flight
76  // are with respect to (0,0,0)
77  avg_time += (rh_time - 1.0);
78  time_norm += 1.0;
79  }
80  if (rh_energy > max_e) {
81  max_e = rh_energy;
82  max_e_layer = rhf.recHitRef()->layer();
83  }
84  if (refhit->detId() == cluster.seed())
85  refseed = refhit;
86  const double rh_fraction = rhf.fraction();
87  rh_energy = refhit->energy() * rh_fraction;
88  if (edm::isNotFinite(rh_energy)) {
89  //temporarily changed exception to warning
90  // throw cms::Exception("PFClusterAlgo")
91  edm::LogWarning("PFClusterAlgo") << "rechit " << refhit->detId() << " has a NaN energy... "
92  << "The input of the particle flow clustering seems to be corrupted.";
93  continue;
94  }
95  pcavars[0] = refhit->position().x();
96  pcavars[1] = refhit->position().y();
97  pcavars[2] = refhit->position().z();
98  int nhit = int(rh_energy * 100); // put rec_hit energy in units of 10 MeV
99 
100  for (int i = 0; i < nhit; ++i) {
101  pca_->AddRow(pcavars);
102  }
103  }
104  cluster.setEnergy(cl_energy);
105  cluster.setLayer(max_e_layer);
106  // calculate the position
107 
108  pca_->MakePrincipals();
109  const TVectorD& means = *(pca_->GetMeanValues());
110  const TMatrixD& eigens = *(pca_->GetEigenVectors());
111 
112  math::XYZPoint barycenter(means[0], means[1], means[2]);
113  math::XYZVector axis(eigens(0, 0), eigens(1, 0), eigens(2, 0));
114 
115  if (time_norm > 0.0) {
116  avg_time = avg_time / time_norm;
117  } else {
119  }
120 
121  if (axis.z() * barycenter.z() < 0.0) {
122  axis = math::XYZVector(-eigens(0, 0), -eigens(1, 0), -eigens(2, 0));
123  }
124 
125  if (updateTiming_) {
126  cluster.setTime(avg_time);
127  }
128  cluster.setPosition(barycenter);
129  cluster.calculatePositionREP();
130 }
void setLayer(PFLayer::Layer layer)
set layer
Definition: PFCluster.cc:49
constexpr bool isNotFinite(T x)
Definition: isFinite.h:9
void setPosition(const math::XYZPoint &p)
Definition: CaloCluster.h:139
DetId seed() const
return DetId of seed
Definition: CaloCluster.h:218
void setTime(float time, float timeError=0)
Definition: PFCluster.h:84
Fraction of a PFRecHit (rechits can be shared between several PFCluster&#39;s)
const std::vector< reco::PFRecHitFraction > & recHitFractions() const
vector of rechit fractions
Definition: PFCluster.h:65
void setEnergy(double energy)
Definition: CaloCluster.h:135
void calculatePositionREP()
computes posrep_ once and for all
Definition: PFCluster.h:95
Layer
layer definition
Definition: PFLayer.h:29
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:31
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
std::unique_ptr< TPrincipal > pca_
Log< level::Warning, false > LogWarning

◆ calculateAndSetPositions()

void Cluster3DPCACalculator::calculateAndSetPositions ( reco::PFClusterCollection clusters,
const HcalPFCuts cuts 
)
overridevirtual

Implements PFCPositionCalculatorBase.

Definition at line 49 of file Cluster3DPCACalculator.cc.

References calculateAndSetPositionActual(), bsc_activity_cfg::clusters, and pca_.

49  {
50  for (reco::PFCluster& cluster : clusters) {
51  pca_ = std::make_unique<TPrincipal>(3, "D");
53  }
54 }
Particle flow cluster, see clustering algorithm in PFClusterAlgo.
Definition: PFCluster.h:42
void calculateAndSetPositionActual(reco::PFCluster &)
std::unique_ptr< TPrincipal > pca_

◆ operator=()

Cluster3DPCACalculator& Cluster3DPCACalculator::operator= ( const Cluster3DPCACalculator )
delete

◆ showerParameters()

void Cluster3DPCACalculator::showerParameters ( const reco::PFCluster ,
math::XYZPoint ,
math::XYZVector  
)
private

Member Data Documentation

◆ pca_

std::unique_ptr<TPrincipal> Cluster3DPCACalculator::pca_
private

◆ updateTiming_

const bool Cluster3DPCACalculator::updateTiming_
private

Definition at line 34 of file Cluster3DPCACalculator.cc.

Referenced by calculateAndSetPositionActual().