CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Types | Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | Friends
PFHcalSuperClusterAlgo Class Reference

Algorithm for particle flow superclustering. More...

#include <PFHcalSuperClusterAlgo.h>

Public Types

typedef edm::Handle
< reco::PFClusterCollection > 		PFClusterHandle
 
typedef edm::Ptr< reco::PFCluster > PFClusterPtr
 

Public Member Functions

std::pair< double, double > calculatePosition (const reco::PFCluster &cluster)
 recalculate eta-phi position of clusters More...
 
std::pair< double, double > calculateWidths (const reco::PFCluster &cluster)
 calculate eta-phi widths of clusters More...
 
const edm::PtrVector
< reco::PFCluster > & 		clusters () const
 
std::auto_ptr< std::vector
< reco::PFCluster > > & 		clusters ()
 setters ----------------------------------------------——— More...
 
void doClustering (const reco::PFClusterCollection &clusters, const reco::PFClusterCollection &clustersHO)
 perform clustering More...
 
void doClustering (const PFClusterHandle &clustersHandle, const PFClusterHandle &clustersHOHandle)
 perform clustering in full framework More...
 
void enableDebugging (bool debug)
 enable/disable debugging More...
 
 PFHcalSuperClusterAlgo ()
 constructor More...
 
std::auto_ptr< std::vector
< reco::PFSuperCluster > > & 		superClusters ()
 
void write ()
 
virtual ~PFHcalSuperClusterAlgo ()
 destructor More...
 

Private Member Functions

void doClusteringWorker (const reco::PFClusterCollection &clusters, const reco::PFClusterCollection &clustersHO)
 perform clustering More...
 

Private Attributes

edm::PtrVector< reco::PFCluster > clusters_
 
PFClusterHandle clustersHandle_
 
PFClusterHandle clustersHOHandle_
 
bool debug_
 debugging on/off More...
 
std::auto_ptr< std::vector
< reco::PFCluster > > 		pfClusters_
 particle flow clusters More...
 
std::auto_ptr< std::vector
< reco::PFSuperCluster > > 		pfSuperClusters_
 particle flow superclusters More...
 
reco::PFSuperCluster SuperCluster_
 

Static Private Attributes

static unsigned prodNum_ = 1
 product number More...
 

Friends

std::ostream & operator<< (std::ostream &out, const PFHcalSuperClusterAlgo &algo)
 

Detailed Description

Algorithm for particle flow superclustering.

This class takes as an input a map of pointers to PFCluster's, and creates PFSuperCluster's from these clusters.

Author
Chris Tully
Date
July 2012

Definition at line 27 of file PFHcalSuperClusterAlgo.h.

Member Typedef Documentation

typedef edm::Handle< reco::PFClusterCollection > PFHcalSuperClusterAlgo::PFClusterHandle

Definition at line 41 of file PFHcalSuperClusterAlgo.h.

typedef edm::Ptr< reco::PFCluster> PFHcalSuperClusterAlgo::PFClusterPtr

Definition at line 42 of file PFHcalSuperClusterAlgo.h.

Constructor & Destructor Documentation

PFHcalSuperClusterAlgo::PFHcalSuperClusterAlgo ( )

constructor

Definition at line 26 of file PFHcalSuperClusterAlgo.cc.

26  :
27  pfClusters_( new vector<reco::PFCluster> ),
28  pfSuperClusters_( new vector<reco::PFSuperCluster> ),
29  debug_(false)
30 {
31 
32 }
PFHcalSuperClusterAlgo::debug_
bool debug_
debugging on/off
Definition: PFHcalSuperClusterAlgo.h:97
PFHcalSuperClusterAlgo::pfClusters_
std::auto_ptr< std::vector< reco::PFCluster > > pfClusters_
particle flow clusters
Definition: PFHcalSuperClusterAlgo.h:85
PFHcalSuperClusterAlgo::pfSuperClusters_
std::auto_ptr< std::vector< reco::PFSuperCluster > > pfSuperClusters_
particle flow superclusters
Definition: PFHcalSuperClusterAlgo.h:88
virtual PFHcalSuperClusterAlgo::~PFHcalSuperClusterAlgo ( )
inlinevirtual

destructor

Definition at line 35 of file PFHcalSuperClusterAlgo.h.

35 {;}

Member Function Documentation

std::pair< double, double > PFHcalSuperClusterAlgo::calculatePosition ( const reco::PFCluster &  cluster)

recalculate eta-phi position of clusters

Definition at line 60 of file PFHcalSuperClusterAlgo.cc.

References cuy::denominator, reco::PFCluster::energy(), create_public_lumi_plots::log, max(), Geom::pi(), reco::PFCluster::recHitFractions(), Geom::twoPi(), and w().

Referenced by calculateWidths(), and doClusteringWorker().

61 {
62  double numeratorEta = 0.0;
63  double numeratorPhi = 0.0;
64  double denominator = 0.0;
65  double posEta = 0.0;
66  double posPhi = 0.0;
67  double w0_ = 4.2;
68  const double clusterEnergy = cluster.energy();
69  if (cluster.energy()>0.0 && cluster.recHitFractions().size()>0) {
70  const std::vector <reco::PFRecHitFraction >& pfhitsandfracs = cluster.recHitFractions();
71  for (std::vector<reco::PFRecHitFraction>::const_iterator it = pfhitsandfracs.begin(); it != pfhitsandfracs.end(); ++it) {
72  const reco::PFRecHitRef rechit = it->recHitRef();
73  double hitEta = rechit->positionREP().Eta();
74  double hitPhi = rechit->positionREP().Phi();
75  //Change into the -pi to +pi angular range
76  while (hitPhi > +Geom::pi()) { hitPhi -= Geom::twoPi(); }
77  while (hitPhi < -Geom::pi()) { hitPhi += Geom::twoPi(); }
78  double hitEnergy = rechit->energy();
79  const double w = std::max(0.0, w0_ + log(hitEnergy / clusterEnergy));
80  denominator += w;
81  numeratorEta += w*hitEta;
82  numeratorPhi += w*hitPhi;
83  }
84  posEta = numeratorEta/denominator;
85  posPhi = numeratorPhi/denominator;
86  }
87 
88  pair<double, double> posEtaPhi(posEta,posPhi);
89 
90  return posEtaPhi;
91 }
create_public_lumi_plots.log
log
Definition: create_public_lumi_plots.py:1119
cuy.denominator
list denominator
Definition: cuy.py:484
max
const T & max(const T &a, const T &b)
Definition: MaterialBudgetTrackerHistos.cc:4
reco::PFCluster::energy
double energy() const
cluster energy
Definition: PFCluster.h:73
Geom::pi
double pi()
Definition: Pi.h:31
Geom::twoPi
double twoPi()
Definition: Pi.h:32
reco::PFCluster::recHitFractions
const std::vector< reco::PFRecHitFraction > & recHitFractions() const
vector of rechit fractions
Definition: PFCluster.h:63
w
T w() const
Definition: extBasic3DVector.h:245
std::pair< double, double > PFHcalSuperClusterAlgo::calculateWidths ( const reco::PFCluster &  cluster)

calculate eta-phi widths of clusters

Definition at line 94 of file PFHcalSuperClusterAlgo.cc.

References abs, calculatePosition(), cuy::denominator, dPhi(), reco::PFCluster::energy(), create_public_lumi_plots::log, max(), Geom::pi(), reco::PFCluster::recHitFractions(), mathSSE::sqrt(), Geom::twoPi(), and w().

Referenced by doClusteringWorker().

95 {
96  double numeratorEtaEta = 0;
97  // double numeratorEtaPhi = 0;
98  double numeratorPhiPhi = 0;
99  double denominator = 0;
100  double widthEta = 0.0;
101  double widthPhi = 0.0;
102 
103  double w0_ = 4.2;
104 
105  const double clusterEta = calculatePosition(cluster).first;
106  const double clusterPhi = calculatePosition(cluster).second;
107  const double clusterEnergy = cluster.energy();
108  if(cluster.energy()>0.0 && cluster.recHitFractions().size()>0) {
109  double hitEta, hitPhi, hitEnergy, dEta, dPhi;
110  const std::vector< reco::PFRecHitFraction >& pfhitsandfracs = cluster.recHitFractions();
111  for (std::vector<reco::PFRecHitFraction>::const_iterator it = pfhitsandfracs.begin(); it != pfhitsandfracs.end(); ++it) {
112  const reco::PFRecHitRef rechit = it->recHitRef();
113  // rechit->calculatePositionREP();
114  hitEta = rechit->positionREP().Eta();
115  hitPhi = rechit->positionREP().Phi();
116  hitEnergy = rechit->energy();
117  dEta = hitEta - clusterEta;
118  dPhi = hitPhi - clusterPhi;
119 
120  while (dPhi > +Geom::pi()) { dPhi -= Geom::twoPi(); }
121  while (dPhi < -Geom::pi()) { dPhi += Geom::twoPi(); }
122 
123 
124  const double w = std::max(0.0, w0_ + log(hitEnergy / clusterEnergy));
125 
126  denominator += w;
127  numeratorEtaEta += w * dEta * dEta;
128  // numeratorEtaPhi += w * dEta * dPhi;
129  numeratorPhiPhi += w * dPhi * dPhi;
130  }
131 
132  double covEtaEta = numeratorEtaEta / denominator;
133  // double covEtaPhi_ = numeratorEtaPhi / denominator;
134  double covPhiPhi = numeratorPhiPhi / denominator;
135 
136  widthEta = sqrt(abs(covEtaEta));
137  widthPhi = sqrt(abs(covPhiPhi));
138  }
139  pair<double, double> widthEtaPhi(widthEta,widthPhi);
140  return widthEtaPhi;
141 }
create_public_lumi_plots.log
log
Definition: create_public_lumi_plots.py:1119
abs
#define abs(x)
Definition: mlp_lapack.h:159
cuy.denominator
list denominator
Definition: cuy.py:484
dPhi
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
max
const T & max(const T &a, const T &b)
Definition: MaterialBudgetTrackerHistos.cc:4
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
reco::PFCluster::energy
double energy() const
cluster energy
Definition: PFCluster.h:73
PFHcalSuperClusterAlgo::calculatePosition
std::pair< double, double > calculatePosition(const reco::PFCluster &cluster)
recalculate eta-phi position of clusters
Definition: PFHcalSuperClusterAlgo.cc:60
Geom::pi
double pi()
Definition: Pi.h:31
Geom::twoPi
double twoPi()
Definition: Pi.h:32
reco::PFCluster::recHitFractions
const std::vector< reco::PFRecHitFraction > & recHitFractions() const
vector of rechit fractions
Definition: PFCluster.h:63
w
T w() const
Definition: extBasic3DVector.h:245
const edm::PtrVector<reco::PFCluster>& PFHcalSuperClusterAlgo::clusters ( ) const
inline

Definition at line 43 of file PFHcalSuperClusterAlgo.h.

References clusters_.

Referenced by doClustering().

44  {return clusters_; }
PFHcalSuperClusterAlgo::clusters_
edm::PtrVector< reco::PFCluster > clusters_
Definition: PFHcalSuperClusterAlgo.h:93
std::auto_ptr< std::vector< reco::PFCluster > >& PFHcalSuperClusterAlgo::clusters ( )
inline

setters ----------------------------------------------———

getters -------------------------------------------------——

Returns
particle flow clusters

Definition at line 69 of file PFHcalSuperClusterAlgo.h.

References pfClusters_.

70  {return pfClusters_;}
PFHcalSuperClusterAlgo::pfClusters_
std::auto_ptr< std::vector< reco::PFCluster > > pfClusters_
particle flow clusters
Definition: PFHcalSuperClusterAlgo.h:85
void PFHcalSuperClusterAlgo::doClustering ( const reco::PFClusterCollection &  clusters,
const reco::PFClusterCollection &  clustersHO 
)

perform clustering

Definition at line 49 of file PFHcalSuperClusterAlgo.cc.

References edm::HandleBase::clear(), clustersHandle_, clustersHOHandle_, and doClusteringWorker().

49  {
50  // using clusters without a Handle, clear to avoid a stale member
51  clustersHandle_.clear();
52  clustersHOHandle_.clear();
53 
54  // perform clustering
55  doClusteringWorker( clusters, clustersHO );
56 }
PFHcalSuperClusterAlgo::doClusteringWorker
void doClusteringWorker(const reco::PFClusterCollection &clusters, const reco::PFClusterCollection &clustersHO)
perform clustering
Definition: PFHcalSuperClusterAlgo.cc:143
PFHcalSuperClusterAlgo::clusters
const edm::PtrVector< reco::PFCluster > & clusters() const
Definition: PFHcalSuperClusterAlgo.h:43
void PFHcalSuperClusterAlgo::doClustering ( const PFClusterHandle &  clustersHandle,
const PFClusterHandle &  clustersHOHandle 
)

perform clustering in full framework

Definition at line 37 of file PFHcalSuperClusterAlgo.cc.

References clusters(), clustersHandle_, clustersHOHandle_, and doClusteringWorker().

37  {
38  const reco::PFClusterCollection& clusters = * clustersHandle;
39  const reco::PFClusterCollection& clustersHO = * clustersHOHandle;
40 
41  // cache the Handle to the clusters
42  clustersHandle_ = clustersHandle;
43  clustersHOHandle_ = clustersHOHandle;
44 
45  // perform clustering
46  doClusteringWorker( clusters, clustersHO );
47 }
PFHcalSuperClusterAlgo::doClusteringWorker
void doClusteringWorker(const reco::PFClusterCollection &clusters, const reco::PFClusterCollection &clustersHO)
perform clustering
Definition: PFHcalSuperClusterAlgo.cc:143
reco::PFClusterCollection
std::vector< PFCluster > PFClusterCollection
collection of PFCluster objects
Definition: PFClusterFwd.h:9
PFHcalSuperClusterAlgo::clusters
const edm::PtrVector< reco::PFCluster > & clusters() const
Definition: PFHcalSuperClusterAlgo.h:43
void PFHcalSuperClusterAlgo::doClusteringWorker ( const reco::PFClusterCollection &  clusters,
const reco::PFClusterCollection &  clustersHO 
)
private

perform clustering

Definition at line 143 of file PFHcalSuperClusterAlgo.cc.

References abs, calculatePosition(), calculateWidths(), edm::PtrVectorBase::clear(), clusters_, clustersHandle_, gather_cfg::cout, SiPixelRawToDigiRegional_cfi::deltaPhi, deltaR(), dPhi(), PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, first, PFLayer::HCAL_BARREL1, PFLayer::HCAL_BARREL2, PFLayer::HCAL_ENDCAP, init, reco::PFHcalSuperClusterInit::initialize(), pfClusters_, pfSuperClusters_, edm::PtrVector< T >::push_back(), edm::PtrVectorBase::size(), mathSSE::sqrt(), w2, w3, and w4.

Referenced by doClustering().

143  {
144 
145  double dRcut=0.17;
146  double dEtacut = 0.0;
147  double dPhicut = 0.0;
148  double etaScale = 1.0;
149  double phiScale = 0.5;
150  // double dRcut=0.30;
151 
152  if ( pfClusters_.get() )
153  pfClusters_->clear();
154  else
155  pfClusters_.reset( new std::vector<reco::PFCluster> );
156 
157  if ( pfSuperClusters_.get() )
158  pfSuperClusters_->clear();
159  else
160  pfSuperClusters_.reset( new std::vector<reco::PFSuperCluster> );
161 
162  // compute cluster depth index
163  std::vector< unsigned > clusterdepth(clusters.size());
164  // cout << " clusters: " << clusters.size() <<endl;
165  int mclustersHB=0;
166  int mclustersHE=0;
167  for (unsigned short ic=0; ic<clusters.size();++ic) {
168  if( clusters[ic].layer() == PFLayer::HCAL_BARREL1) mclustersHB++;
169  if( clusters[ic].layer() == PFLayer::HCAL_ENDCAP) mclustersHE++;
170  }
171  for (unsigned short ic=0; ic<clusters.size();++ic)
172  {
173  if( clusters[ic].layer() == PFLayer::HCAL_BARREL1
174  || clusters[ic].layer() == PFLayer::HCAL_ENDCAP ) { //Hcal case
175 
176  const std::vector< std::pair<DetId, float> > & hitsandfracs =
177  clusters[ic].hitsAndFractions();
178  unsigned clusterdepthfirst=0;
179  for(unsigned ihandf=0; ihandf<hitsandfracs.size(); ihandf++) {
180  unsigned depth = ((HcalDetId)hitsandfracs[ihandf].first).depth();
181  // cout << " depth parameter from clustering: " << depth <<endl;
182  clusterdepth[ic] = depth;
183  if( ihandf>0 && depth!=clusterdepthfirst) cout << " Problem with cluster depth: " << depth << " not equal to " << clusterdepthfirst <<endl;
184  else if( ihandf==0 ) clusterdepthfirst = depth;
185  }
186  // delete hitsandfracs;
187  }
188  }
189  std::vector< unsigned > clusterdepthHO(clustersHO.size());
190  // cout << " HO clusters: " << clustersHO.size() <<endl;
191  double hcaleta1=0.0;
192  double hcalphi1=0.0;
193  double hcaleta2=0.0;
194  double hcalphi2=0.0;
195  double dR = 0.0;
196  double dEta = 0.0;
197  double dPhi = 0.0;
198  for (unsigned short ic=0; ic<clustersHO.size();++ic)
199  {
200  if( clustersHO[ic].layer() == PFLayer::HCAL_BARREL2) { //HO case
201 
202  const std::vector< std::pair<DetId, float> > & hitsandfracs =
203  clustersHO[ic].hitsAndFractions();
204  unsigned clusterdepthfirst=0;
205  for(unsigned ihandf=0; ihandf<hitsandfracs.size(); ihandf++) {
206  unsigned depth = ((HcalDetId)hitsandfracs[ihandf].first).depth();
207  // cout << " depth parameter from HO clustering: " << depth <<endl;
208  clusterdepthHO[ic] = depth;
209  if( ihandf>0 && depth!=clusterdepthfirst) cout << " Problem with HO cluster depth: " << depth << " not equal to " << clusterdepthfirst <<endl;
210  else if( ihandf==0 ) clusterdepthfirst = depth;
211  }
212  // delete hitsandfracs;
213  }
214  }
215 
216  std::vector< unsigned > imerge(clusters.size());
217  std::vector< unsigned > imergeHO(clustersHO.size());
218  std::vector< bool > lmerge(clusters.size());
219  std::vector< bool > lmergeHO(clustersHO.size());
220 
221  hcaleta1=0.0;
222  hcalphi1=0.0;
223  hcaleta2=0.0;
224  hcalphi2=0.0;
225  dR = 0.0;
226  dEta = 0.0;
227  dPhi = 0.0;
228 
229  // cout << " setting up cluster merging indices "<<endl;
230  for (unsigned short ic1=0; ic1<clusters.size();++ic1) {
231  lmerge[ic1]=false;
232  }
233  for (unsigned short ic1=0; ic1<clustersHO.size();++ic1) {
234  lmergeHO[ic1]=false;
235  }
236  for (unsigned short ic1=0; ic1<clusters.size();++ic1) {
237  if( clusterdepth[ic1]==1 ){
238  hcaleta1 = calculatePosition(clusters[ic1]).first;
239  hcalphi1 = calculatePosition(clusters[ic1]).second;
240  for (unsigned short ic2=0; ic2<clusters.size();++ic2) {
241  hcaleta2 = calculatePosition(clusters[ic2]).first;
242  hcalphi2 = calculatePosition(clusters[ic2]).second;
243  dR = deltaR( hcaleta1, hcalphi1, hcaleta2, hcalphi2 );
244  dEta = abs(hcaleta1 - hcaleta2);
245  dPhi = abs(deltaPhi(hcalphi1, hcalphi2));
246  double w1 = calculateWidths(clusters[ic1]).first;
247  if (w1 == 0) w1 = 0.087;
248  double w2 = calculateWidths(clusters[ic2]).first;
249  if (w2 == 0) w2 = 0.087;
250  double w3 = calculateWidths(clusters[ic1]).second;
251  if (w3 < 0.087) w3 = 0.087;
252  double w4 = calculateWidths(clusters[ic2]).second;
253  if (w4 < 0.087) w4 = 0.087;
254  double etawidth = sqrt(w1*w1 + w2*w2);
255  double phiwidth = sqrt(w3*w3 + w4*w4);
256  dEtacut = etaScale*etawidth;
257  dPhicut = phiScale*phiwidth;
258  if( clusterdepth[ic2]==2 ){
259  // cout << " depth 1-2 dR = " << dR <<endl;
260  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
261  imerge[ic2]=ic1;
262  lmerge[ic2]=true;
263  }
264  } else if( clusterdepth[ic2]==3 ){
265  // cout << " depth 1-3 dR = " << dR <<endl;
266  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut) ) {
267  imerge[ic2]=ic1;
268  lmerge[ic2]=true;
269  }
270  }
271  }
272  } else if( clusterdepth[ic1]==2 ){
273  hcaleta1 = calculatePosition(clusters[ic1]).first;
274  hcalphi1 = calculatePosition(clusters[ic1]).second;
275  for (unsigned short ic2=0; ic2<clusters.size();++ic2) {
276  hcaleta2 = calculatePosition(clusters[ic2]).first;
277  hcalphi2 = calculatePosition(clusters[ic2]).second;
278  dEta = abs(hcaleta1 - hcaleta2);
279  dPhi = abs(deltaPhi(hcalphi1, hcalphi2));
280  dR = deltaR( hcaleta1, hcalphi1, hcaleta2, hcalphi2 );
281  double w1 = calculateWidths(clusters[ic1]).first;
282  if (w1 == 0) w1 = 0.087;
283  double w2 = calculateWidths(clusters[ic2]).first;
284  if (w2 == 0) w2 = 0.087;
285  double w3 = calculateWidths(clusters[ic1]).second;
286  if (w3 < 0.087) w3 = 0.087;
287  double w4 = calculateWidths(clusters[ic2]).second;
288  if (w4 < 0.087) w4 = 0.087;
289  double etawidth = sqrt(w1*w1+w2*w2);
290  double phiwidth = sqrt(w3*w3+w4*w4);
291  dEtacut = etaScale*etawidth;
292  dPhicut = phiScale*phiwidth;
293  if( clusterdepth[ic2]==3 ){
294  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
295  imerge[ic2]=ic1;
296  lmerge[ic2]=true;
297  }
298  } else if( clusterdepth[ic2]==4 ){
299  // cout << " depth 2-4 dR = " << dR <<endl;
300  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
301  imerge[ic2]=ic1;
302  lmerge[ic2]=true;
303  }
304  }
305  }
306  } else if( clusterdepth[ic1]==3 ){
307  hcaleta1 = calculatePosition(clusters[ic1]).first;
308  hcalphi1 = calculatePosition(clusters[ic1]).second;
309  for (unsigned short ic2=0; ic2<clusters.size();++ic2) {
310  hcaleta2 = calculatePosition(clusters[ic2]).first;
311  hcalphi2 = calculatePosition(clusters[ic2]).second;
312  dEta = abs(hcaleta1 - hcaleta2);
313  dPhi = abs(deltaPhi(hcalphi1, hcalphi2));
314  dR = deltaR( hcaleta1, hcalphi1, hcaleta2, hcalphi2 );
315  double w1 = calculateWidths(clusters[ic1]).first;
316  if (w1 == 0) w1 = 0.087;
317  double w2 = calculateWidths(clusters[ic2]).first;
318  if (w2 == 0) w2 = 0.087;
319  double w3 = calculateWidths(clusters[ic1]).second;
320  if (w3 < 0.087) w3 = 0.087;
321  double w4 = calculateWidths(clusters[ic2]).second;
322  if (w4 < 0.087) w4 = 0.087;
323  double etawidth = sqrt(w1*w1+w2*w2);
324  double phiwidth = sqrt(w3*w3+w4*w4);
325  dEtacut = etaScale*etawidth;
326  dPhicut = phiScale*phiwidth;
327  if( clusterdepth[ic2]==4 ){
328  // cout << " depth 3-4 dR = " << dR <<endl;
329  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
330  imerge[ic2]=ic1;
331  lmerge[ic2]=true;
332  }
333  } else if( clusterdepth[ic2]==5 ){
334  // cout << " depth 3-5 dR = " << dR <<endl;
335  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
336  imerge[ic2]=ic1;
337  lmerge[ic2]=true;
338  }
339  }
340  }
341  for (unsigned short ic2=0; ic2<clustersHO.size();++ic2) {
342  hcaleta2 = calculatePosition(clustersHO[ic2]).first;
343  hcalphi2 = calculatePosition(clustersHO[ic2]).second;
344  dEta = abs(hcaleta1 - hcaleta2);
345  dPhi = abs(deltaPhi(hcalphi1, hcalphi2));
346  dR = deltaR( hcaleta1, hcalphi1, hcaleta2, hcalphi2 );
347  double w1 = calculateWidths(clusters[ic1]).first;
348  if (w1 == 0) w1 = 0.087;
349  double w2 = calculateWidths(clustersHO[ic2]).first;
350  if (w2 == 0) w2 = 0.087;
351  double w3 = calculateWidths(clusters[ic1]).second;
352  if (w3 < 0.087) w3 = 0.087;
353  double w4 = calculateWidths(clustersHO[ic2]).second;
354  if (w4 < 0.087) w4 = 0.087;
355  double etawidth = sqrt(w1*w1+w2*w2);
356  double phiwidth = sqrt(w3*w3+w4*w4);
357  dEtacut = etaScale*etawidth;
358  dPhicut = phiScale*phiwidth;
359  if( clusterdepthHO[ic2]==5 ){
360  // cout << " depth 3-HO dR = " << dR <<endl;
361  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
362  imergeHO[ic2]=ic1;
363  lmergeHO[ic2]=true;
364  }
365  }
366  }
367  } else if( clusterdepth[ic1]==4 ){
368  hcaleta1 = calculatePosition(clusters[ic1]).first;
369  hcalphi1 = calculatePosition(clusters[ic1]).second;
370  for (unsigned short ic2=0; ic2<clusters.size();++ic2) {
371  hcaleta2 = calculatePosition(clusters[ic2]).first;
372  hcalphi2 = calculatePosition(clusters[ic2]).second;
373  dEta = abs(hcaleta1 - hcaleta2);
374  dPhi = abs(deltaPhi(hcalphi1, hcalphi2));
375  dR = deltaR( hcaleta1, hcalphi1, hcaleta2, hcalphi2 );
376  double w1 = calculateWidths(clusters[ic1]).first;
377  if (w1 < 0.087) w1 = 0.087;
378  double w2 = calculateWidths(clusters[ic2]).first;
379  if (w2 < 0.087) w2 = 0.087;
380  double w3 = calculateWidths(clusters[ic1]).second;
381  if (w3 < 0.087) w3 = 0.087;
382  double w4 = calculateWidths(clusters[ic2]).second;
383  if (w4 < 0.087) w4 = 0.087;
384  double etawidth = sqrt(w1*w1+w2*w2);
385  double phiwidth = sqrt(w3*w3+w4*w4);
386  dEtacut = etaScale*etawidth;
387  dPhicut = phiScale*phiwidth;
388  if( clusterdepth[ic2]==5 ){
389  // cout << " depth 4-5 dR = " << dR <<endl;
390  if((dR<dRcut) || (dEta<dEtacut && dPhi<dPhicut)) {
391  imerge[ic2]=ic1;
392  lmerge[ic2]=true;
393  }
394  }
395  }
396  } else if( clusterdepth[ic1]==5 ){
397  hcaleta1 = calculatePosition(clusters[ic1]).first;
398  hcalphi1 = calculatePosition(clusters[ic1]).second;
399  }
400  }
401 
402  // start a supercluster with a depth=1 cluster, then loop on all other
403  // clusters to add to cluster list, then for each cluster to add, loop
404  // on remaining clusters to check 2nd level of addition, repeat for all layers
405 
406  // need to add HO cluster logic
407  edm::PtrVector< reco::PFCluster > mergeclusters;
408  std::vector< bool > lmergeclusters(clusters.size());
409  for (unsigned short id=0; id<4;++id) {
410  // cout << " merging with starting depth: "<<id<<endl;
411  for (unsigned short ic1=0; ic1<clusters.size();++ic1) {
412  if( clusterdepth[ic1]==(unsigned)(1+id) ){
413  if(!lmerge[ic1]) {
414  for (unsigned short ic=0; ic<clusters.size();++ic) {
415  lmergeclusters[ic]=false;
416  }
417  // mergeclusters.push_back(clusters[ic1]);
418  for (unsigned short ic2=0; ic2<clusters.size();++ic2) {
419  if( clusterdepth[ic2]==(unsigned)(2+id) ){
420  if( imerge[ic2]==ic1 && lmerge[ic2] ) {
421  // mergeclusters.push_back(clusters[ic2]);
422  lmergeclusters[ic2]=true;
423  for (unsigned short ic3=0; ic3<clusters.size();++ic3) {
424  if( clusterdepth[ic3]==(unsigned)(3+id) ){
425  if( imerge[ic3]==ic2 && lmerge[ic3] ) {
426  // mergeclusters.push_back(clusters[ic3]);
427  lmergeclusters[ic3]=true;
428  for (unsigned short ic4=0; ic4<clusters.size();++ic4) {
429  if( clusterdepth[ic4]==(unsigned)(4+id) ){
430  if( imerge[ic4]==ic3 && lmerge[ic4] ) {
431  // mergeclusters.push_back(clusters[ic4]);
432  lmergeclusters[ic4]=true;
433  for (unsigned short ic5=0; ic5<clusters.size();++ic5) {
434  if( clusterdepth[ic5]==(unsigned)(5+id) ){
435  // if( imerge[ic5]==ic4 && lmerge[ic5] ) mergeclusters.push_back(clusters[ic5]);
436  if( imerge[ic5]==ic4 && lmerge[ic5] ) lmergeclusters[ic5]=true;
437  }
438  }
439  } else if( clusterdepth[ic4]==(unsigned)(5+id) ){
440  // if( imerge[ic4]==ic3 && lmerge[ic4] ) mergeclusters.push_back(clusters[ic4]);
441  if( imerge[ic4]==ic3 && lmerge[ic4] ) lmergeclusters[ic4]=true;
442  }
443  }
444  }
445  } else if( clusterdepth[ic3]==(unsigned)(4+id) ){
446  if( imerge[ic3]==ic2 && lmerge[ic3] ) {
447  // mergeclusters.push_back(clusters[ic3]);
448  lmergeclusters[ic3]=true;
449  for (unsigned short ic4=0; ic4<clusters.size();++ic4) {
450  if( clusterdepth[ic4]==(unsigned)(5+id) ){
451  // if( imerge[ic4]==ic3 && lmerge[ic4] ) mergeclusters.push_back(clusters[ic4]);
452  if( imerge[ic4]==ic3 && lmerge[ic4] ) lmergeclusters[ic4]=true;
453  }
454  }
455  }
456  }
457  } if( clusterdepth[ic3]==(unsigned)(5+id) ){
458  // if( imerge[ic3]==ic2 && lmerge[ic3] ) mergeclusters.push_back(clusters[ic3]);
459  if( imerge[ic3]==ic2 && lmerge[ic3] ) lmergeclusters[ic3]=true;
460  }
461  }
462  }
463  } else if( clusterdepth[ic2]==(unsigned)(3+id) ){
464  if( imerge[ic2]==ic1 && lmerge[ic2] ) {
465  // mergeclusters.push_back(clusters[ic2]);
466  lmergeclusters[ic2]=true;
467  for (unsigned short ic3=0; ic3<clusters.size();++ic3) {
468  if( clusterdepth[ic3]==(unsigned)(4+id) ){
469  if( imerge[ic3]==ic2 && lmerge[ic3] ) {
470  // mergeclusters.push_back(clusters[ic3]);
471  lmergeclusters[ic3]=true;
472  for (unsigned short ic4=0; ic4<clusters.size();++ic4) {
473  if( clusterdepth[ic4]==(unsigned)(5+id) ){
474  // if( imerge[ic4]==ic3 && lmerge[ic4] ) mergeclusters.push_back(clusters[ic4]);
475  if( imerge[ic4]==ic3 && lmerge[ic4] ) lmergeclusters[ic4]=true;
476  }
477  }
478  } else if( clusterdepth[ic3]==(unsigned)(5+id) ){
479  // if( imerge[ic3]==ic2 && lmerge[ic3] ) mergeclusters.push_back(clusters[ic3]);
480  if( imerge[ic3]==ic2 && lmerge[ic3] ) lmergeclusters[ic3]=true;
481  }
482  }
483  }
484  }
485  }
486  }
487  mergeclusters.push_back( PFClusterPtr( clustersHandle_, ic1));
488  for (unsigned short ic=0; ic<clusters.size();++ic) {
489  if(lmergeclusters[ic]) {
490  mergeclusters.push_back( PFClusterPtr(clustersHandle_, ic ) );
491  }
492  }
493  if(mergeclusters.size()>0) {
494  // cout << " number of clusters to merge: " <<mergeclusters.size()<<endl;
495  reco::PFSuperCluster ipfsupercluster(mergeclusters);
496  PFHcalSuperClusterInit init;
497  init.initialize( ipfsupercluster, clusters_);
498  pfSuperClusters_->push_back(ipfsupercluster);
499  pfClusters_->push_back((reco::PFCluster)ipfsupercluster);
500  mergeclusters.clear();
501  }
502 
503  }
504  } // end of depth 1+id initiated logic
505  }
506  }
507 
508  /*
509  for (unsigned short ic=0; ic<clusters.size();++ic) {
510  mergeclusters.clear();
511  mergeclusters.push_back(clusters[ic]);
512  reco::PFSuperCluster ipfsupercluster(mergeclusters);
513  pfSuperClusters_->push_back(ipfsupercluster);
514  pfClusters_->push_back((reco::PFCluster)ipfsupercluster);
515  // pfClusters_->push_back(clusters[ic]);
516  }
517  */
518 
519  clusterdepth.clear();
520  clusterdepthHO.clear();
521  imerge.clear();
522  imergeHO.clear();
523  lmerge.clear();
524  lmergeHO.clear();
525  mergeclusters.clear();
526 
527 }
edm::PtrVectorBase::size
size_type size() const
Size of the RefVector.
Definition: PtrVectorBase.h:73
w2
common ppss p3p6s2 common epss epspn46 common const1 w2
Definition: inclppp.h:1
reco::PFCluster
Particle flow cluster, see clustering algorithm in PFClusterAlgo.
Definition: PFCluster.h:42
edm::PtrVector::push_back
void push_back(Ptr< T > const &iPtr)
Definition: PtrVector.h:137
init
int init
Definition: HydjetWrapper.h:63
reco::PFSuperCluster
Particle flow cluster, see clustering algorithm in PFSuperClusterAlgo.
Definition: PFSuperCluster.h:26
abs
#define abs(x)
Definition: mlp_lapack.h:159
PFHcalSuperClusterAlgo::clusters_
edm::PtrVector< reco::PFCluster > clusters_
Definition: PFHcalSuperClusterAlgo.h:93
w4
common ppss p3p6s2 common epss epspn46 common const1 w4
Definition: inclppp.h:1
dPhi
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
PFHcalSuperClusterAlgo::PFClusterPtr
edm::Ptr< reco::PFCluster > PFClusterPtr
Definition: PFHcalSuperClusterAlgo.h:42
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi.dR
tuple dR
Definition: PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi.py:19
PFHcalSuperClusterAlgo::calculateWidths
std::pair< double, double > calculateWidths(const reco::PFCluster &cluster)
calculate eta-phi widths of clusters
Definition: PFHcalSuperClusterAlgo.cc:94
first
bool first
Definition: L1TdeRCT.cc:94
reco::PFHcalSuperClusterInit::initialize
void initialize(reco::PFSuperCluster &supercluster, edm::PtrVector< reco::PFCluster > const &clusters)
Definition: PFHcalSuperClusterInit.cc:10
deltaR
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
PFHcalSuperClusterAlgo::pfClusters_
std::auto_ptr< std::vector< reco::PFCluster > > pfClusters_
particle flow clusters
Definition: PFHcalSuperClusterAlgo.h:85
SiPixelRawToDigiRegional_cfi.deltaPhi
tuple deltaPhi
Definition: SiPixelRawToDigiRegional_cfi.py:9
PFHcalSuperClusterAlgo::calculatePosition
std::pair< double, double > calculatePosition(const reco::PFCluster &cluster)
recalculate eta-phi position of clusters
Definition: PFHcalSuperClusterAlgo.cc:60
w3
common ppss p3p6s2 common epss epspn46 common const1 w3
Definition: inclppp.h:1
edm::PtrVectorBase::clear
void clear()
Clear the PtrVector.
Definition: PtrVectorBase.h:79
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
PFHcalSuperClusterAlgo::clusters
const edm::PtrVector< reco::PFCluster > & clusters() const
Definition: PFHcalSuperClusterAlgo.h:43
PFHcalSuperClusterAlgo::pfSuperClusters_
std::auto_ptr< std::vector< reco::PFSuperCluster > > pfSuperClusters_
particle flow superclusters
Definition: PFHcalSuperClusterAlgo.h:88
void PFHcalSuperClusterAlgo::enableDebugging ( bool  debug)
inline

enable/disable debugging

Definition at line 38 of file PFHcalSuperClusterAlgo.h.

References debug, and debug_.

38 { debug_ = debug;}
PFHcalSuperClusterAlgo::debug_
bool debug_
debugging on/off
Definition: PFHcalSuperClusterAlgo.h:97
debug
#define debug
Definition: MEtoEDMFormat.h:34
std::auto_ptr< std::vector< reco::PFSuperCluster > >& PFHcalSuperClusterAlgo::superClusters ( )
inline
Returns
particle flow superclusters

Definition at line 73 of file PFHcalSuperClusterAlgo.h.

References pfSuperClusters_.

74  {return pfSuperClusters_;}
PFHcalSuperClusterAlgo::pfSuperClusters_
std::auto_ptr< std::vector< reco::PFSuperCluster > > pfSuperClusters_
particle flow superclusters
Definition: PFHcalSuperClusterAlgo.h:88
void PFHcalSuperClusterAlgo::write ( void  )

Definition at line 35 of file PFHcalSuperClusterAlgo.cc.

35  {
36 }

Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  out,
const PFHcalSuperClusterAlgo &  algo 
)
friend

Member Data Documentation

edm::PtrVector< reco::PFCluster > PFHcalSuperClusterAlgo::clusters_
private

Definition at line 93 of file PFHcalSuperClusterAlgo.h.

Referenced by clusters(), and doClusteringWorker().

PFClusterHandle PFHcalSuperClusterAlgo::clustersHandle_
private

Definition at line 91 of file PFHcalSuperClusterAlgo.h.

Referenced by doClustering(), and doClusteringWorker().

PFClusterHandle PFHcalSuperClusterAlgo::clustersHOHandle_
private

Definition at line 92 of file PFHcalSuperClusterAlgo.h.

Referenced by doClustering().

bool PFHcalSuperClusterAlgo::debug_
private

debugging on/off

Definition at line 97 of file PFHcalSuperClusterAlgo.h.

Referenced by enableDebugging().

std::auto_ptr< std::vector<reco::PFCluster> > PFHcalSuperClusterAlgo::pfClusters_
private

particle flow clusters

Definition at line 85 of file PFHcalSuperClusterAlgo.h.

Referenced by clusters(), doClusteringWorker(), and operator<<().

std::auto_ptr< std::vector<reco::PFSuperCluster> > PFHcalSuperClusterAlgo::pfSuperClusters_
private

particle flow superclusters

Definition at line 88 of file PFHcalSuperClusterAlgo.h.

Referenced by doClusteringWorker(), operator<<(), and superClusters().

unsigned PFHcalSuperClusterAlgo::prodNum_ = 1
staticprivate

product number

Definition at line 100 of file PFHcalSuperClusterAlgo.h.

reco::PFSuperCluster PFHcalSuperClusterAlgo::SuperCluster_
private

Definition at line 90 of file PFHcalSuperClusterAlgo.h.

Generated for CMSSW Reference Manual by doxygen 1.8.5