CMS 3D CMS Logo
List of all members | Public Member Functions | Private Member Functions | Private Attributes
HemisphereAlgo Class Reference

#include <HemisphereAlgo.h>

Public Member Functions

std::vector< float > getAxis1 ()
 
std::vector< float > getAxis2 ()
 
std::vector< int > getGrouping ()
 
 HemisphereAlgo (const std::vector< reco::CandidatePtr > &componentRefs_, const int seed_method=0, const int hemisphere_association_method=0)
 
void SetMethod (int seed_method, int hemisphere_association_method)
 
void SetNoSeed (int object_number)
 
 ~HemisphereAlgo ()
 

Private Member Functions

int reconstruct ()
 

Private Attributes

std::vector< float > Axis1
 
std::vector< float > Axis2
 
int hemi_meth
 
std::vector< reco::CandidatePtr > Object
 
std::vector< int > Object_Group
 
std::vector< int > Object_Noseed
 
int seed_meth
 
int status
 

Detailed Description

Definition at line 23 of file HemisphereAlgo.h.

Constructor & Destructor Documentation

◆ HemisphereAlgo()

HemisphereAlgo::HemisphereAlgo ( const std::vector< reco::CandidatePtr > &  componentRefs_,
const int  seed_method = 0,
const int  hemisphere_association_method = 0 
)

Definition at line 11 of file HemisphereAlgo.cc.

14  : Object(componentPtr_),
15  Object_Group(),
16  Axis1(),
17  Axis2(),
18  seed_meth(seed_method),
19  hemi_meth(hemisphere_association_method),
20  status(0) {
21  for (int i = 0; i < (int)Object.size(); i++) {
22  Object_Noseed.push_back(0);
23  }
24 }

References mps_fire::i, createfilelist::int, Object, and Object_Noseed.

◆ ~HemisphereAlgo()

HemisphereAlgo::~HemisphereAlgo ( )
inline

Definition at line 52 of file HemisphereAlgo.h.

52 {};

Member Function Documentation

◆ getAxis1()

vector< float > HemisphereAlgo::getAxis1 ( )

Definition at line 26 of file HemisphereAlgo.cc.

26  {
27  if (status != 1) {
28  this->reconstruct();
29  }
30  return Axis1;
31 }

References Axis1, reconstruct(), and status.

Referenced by PATHemisphereProducer::produce().

◆ getAxis2()

vector< float > HemisphereAlgo::getAxis2 ( )

Definition at line 32 of file HemisphereAlgo.cc.

32  {
33  if (status != 1) {
34  this->reconstruct();
35  }
36  return Axis2;
37 }

References Axis2, reconstruct(), and status.

Referenced by PATHemisphereProducer::produce().

◆ getGrouping()

vector< int > HemisphereAlgo::getGrouping ( )

Definition at line 39 of file HemisphereAlgo.cc.

39  {
40  if (status != 1) {
41  this->reconstruct();
42  }
43  return Object_Group;
44 }

References Object_Group, reconstruct(), and status.

Referenced by PATHemisphereProducer::produce().

◆ reconstruct()

int HemisphereAlgo::reconstruct ( )
private

Definition at line 46 of file HemisphereAlgo.cc.

46  {
47  int vsize = (int)Object.size();
48 
49  LogDebug("HemisphereAlgo") << " HemisphereAlgo method ";
50 
51  Object_Group.clear();
52  Axis1.clear();
53  Axis2.clear();
54 
55  for (int j = 0; j < vsize; j++) {
56  Object_Group.push_back(0);
57  }
58 
59  for (int j = 0; j < 5; j++) {
60  Axis1.push_back(0);
61  Axis2.push_back(0);
62  }
63 
64  for (int i = 0; i < vsize; i++) {
65  if ((*(Object)[i]).p() > (*Object[i]).energy() + 0.001) {
66  edm::LogWarning("HemisphereAlgo") << "Object " << i << " has E = " << (*Object[i]).energy()
67  << " less than P = " << (*Object[i]).p();
68  }
69  }
70 
71  LogDebug("HemisphereAlgo") << " Seeding method = " << seed_meth;
72  int I_Max = -1;
73  int J_Max = -1;
74 
75  if (seed_meth == 1) {
76  float P_Max = 0.;
77  float DeltaRP_Max = 0.;
78 
79  // take highest momentum object as first axis
80  for (int i = 0; i < vsize; i++) {
81  Object_Group[i] = 0;
82  if (Object_Noseed[i] == 0 && P_Max < (*Object[i]).p()) {
83  P_Max = (*Object[i]).p();
84  I_Max = i;
85  }
86  }
87 
88  Axis1[0] = (*Object[I_Max]).px() / (*Object[I_Max]).p();
89  Axis1[1] = (*Object[I_Max]).py() / (*Object[I_Max]).p();
90  Axis1[2] = (*Object[I_Max]).pz() / (*Object[I_Max]).p();
91  Axis1[3] = (*Object[I_Max]).p();
92  Axis1[4] = (*Object[I_Max]).energy();
93 
94  // take as second axis
95  for (int i = 0; i < vsize; i++) {
96  float DeltaR = deltaR((*Object[i]).eta(), (*Object[i]).phi(), (*Object[I_Max]).eta(), (*Object[I_Max]).phi());
97 
98  if (Object_Noseed[i] == 0 && DeltaR > 0.5) {
99  float DeltaRP = DeltaR * (*Object[i]).p();
100  if (DeltaRP > DeltaRP_Max) {
101  DeltaRP_Max = DeltaRP;
102  J_Max = i;
103  }
104  }
105  }
106 
107  if (J_Max >= 0) {
108  Axis2[0] = (*Object[J_Max]).px() / (*Object[J_Max]).p();
109  Axis2[1] = (*Object[J_Max]).py() / (*Object[J_Max]).p();
110  Axis2[2] = (*Object[J_Max]).pz() / (*Object[J_Max]).p();
111  Axis2[3] = (*Object[J_Max]).p();
112  Axis2[4] = (*Object[J_Max]).energy();
113 
114  } else {
115  return 0;
116  }
117  LogDebug("HemisphereAlgo") << " Axis 1 is Object = " << I_Max;
118  LogDebug("HemisphereAlgo") << " Axis 2 is Object = " << J_Max;
119 
120  } else if ((seed_meth == 2) | (seed_meth == 3)) {
121  float Mass_Max = 0.;
122  float InvariantMass = 0.;
123 
124  // maximize the invariant mass of two objects
125  for (int i = 0; i < vsize; i++) {
126  Object_Group[i] = 0;
127  if (Object_Noseed[i] == 0) {
128  for (int j = i + 1; j < vsize; j++) {
129  if (Object_Noseed[j] == 0) {
130  // either the invariant mass
131  if (seed_meth == 2) {
132  InvariantMass =
133  ((*Object[i]).energy() + (*Object[j]).energy()) * ((*Object[i]).energy() + (*Object[j]).energy()) -
134  ((*Object[i]).px() + (*Object[j]).px()) * ((*Object[i]).px() + (*Object[j]).px()) -
135  ((*Object[i]).py() + (*Object[j]).py()) * ((*Object[i]).py() + (*Object[j]).py()) -
136  ((*Object[i]).pz() + (*Object[j]).pz()) * ((*Object[i]).pz() + (*Object[j]).pz());
137  }
138  // or the transverse mass
139  else if (seed_meth == 3) {
140  float pti = sqrt((*Object[i]).px() * (*Object[i]).px() + (*Object[i]).py() * (*Object[i]).py());
141  float ptj = sqrt((*Object[j]).px() * (*Object[j]).px() + (*Object[j]).py() * (*Object[j]).py());
142  InvariantMass =
143  2. * (pti * ptj - (*Object[i]).px() * (*Object[j]).px() - (*Object[i]).py() * (*Object[j]).py());
144  }
145  if (Mass_Max < InvariantMass) {
146  Mass_Max = InvariantMass;
147  I_Max = i;
148  J_Max = j;
149  }
150  }
151  }
152  }
153  }
154 
155  if (J_Max > 0) {
156  Axis1[0] = (*Object[I_Max]).px() / (*Object[I_Max]).p();
157  Axis1[1] = (*Object[I_Max]).py() / (*Object[I_Max]).p();
158  Axis1[2] = (*Object[I_Max]).pz() / (*Object[I_Max]).p();
159 
160  Axis1[3] = (*Object[I_Max]).p();
161  Axis1[4] = (*Object[I_Max]).energy();
162 
163  Axis2[0] = (*Object[J_Max]).px() / (*Object[J_Max]).p();
164  Axis2[1] = (*Object[J_Max]).py() / (*Object[J_Max]).p();
165  Axis2[2] = (*Object[J_Max]).pz() / (*Object[J_Max]).p();
166 
167  Axis2[3] = (*Object[J_Max]).p();
168  Axis2[4] = (*Object[J_Max]).energy();
169 
170  } else {
171  return 0;
172  }
173  LogDebug("HemisphereAlgo") << " Axis 1 is Object = " << I_Max;
174  LogDebug("HemisphereAlgo") << " Axis 2 is Object = " << J_Max;
175 
176  } else {
177  throw cms::Exception("Configuration") << "Please give a valid seeding method!";
178  }
179 
180  // seeding done
181  // now do the hemisphere association
182 
183  LogDebug("HemisphereAlgo") << " Association method = " << hemi_meth;
184 
185  int numLoop = 0;
186  bool I_Move = true;
187 
188  while (I_Move && (numLoop < 100)) {
189  I_Move = false;
190  numLoop++;
191  LogDebug("HemisphereAlgo") << " Iteration = " << numLoop;
192 
193  float Sum1_Px = 0.;
194  float Sum1_Py = 0.;
195  float Sum1_Pz = 0.;
196  float Sum1_P = 0.;
197  float Sum1_E = 0.;
198  float Sum2_Px = 0.;
199  float Sum2_Py = 0.;
200  float Sum2_Pz = 0.;
201  float Sum2_P = 0.;
202  float Sum2_E = 0.;
203 
204  if (hemi_meth == 1) {
205  for (int i = 0; i < vsize; i++) {
206  float P_Long1 = (*Object[i]).px() * Axis1[0] + (*Object[i]).py() * Axis1[1] + (*Object[i]).pz() * Axis1[2];
207  float P_Long2 = (*Object[i]).px() * Axis2[0] + (*Object[i]).py() * Axis2[1] + (*Object[i]).pz() * Axis2[2];
208  if (P_Long1 >= P_Long2) {
209  if (Object_Group[i] != 1) {
210  I_Move = true;
211  }
212  Object_Group[i] = 1;
213  Sum1_Px += (*Object[i]).px();
214  Sum1_Py += (*Object[i]).py();
215  Sum1_Pz += (*Object[i]).pz();
216  Sum1_P += (*Object[i]).p();
217  Sum1_E += (*Object[i]).energy();
218  } else {
219  if (Object_Group[i] != 2) {
220  I_Move = true;
221  }
222  Object_Group[i] = 2;
223  Sum2_Px += (*Object[i]).px();
224  Sum2_Py += (*Object[i]).py();
225  Sum2_Pz += (*Object[i]).pz();
226  Sum2_P += (*Object[i]).p();
227  Sum2_E += (*Object[i]).energy();
228  }
229  }
230 
231  } else if (hemi_meth == 2 || hemi_meth == 3) {
232  for (int i = 0; i < vsize; i++) {
233  if (i == I_Max) {
234  Object_Group[i] = 1;
235  Sum1_Px += (*Object[i]).px();
236  Sum1_Py += (*Object[i]).py();
237  Sum1_Pz += (*Object[i]).pz();
238  Sum1_P += (*Object[i]).p();
239  Sum1_E += (*Object[i]).energy();
240  } else if (i == J_Max) {
241  Object_Group[i] = 2;
242  Sum2_Px += (*Object[i]).px();
243  Sum2_Py += (*Object[i]).py();
244  Sum2_Pz += (*Object[i]).pz();
245  Sum2_P += (*Object[i]).p();
246  Sum2_E += (*Object[i]).energy();
247  } else {
248  if (!I_Move) {
249  float NewAxis1_Px = Axis1[0] * Axis1[3];
250  float NewAxis1_Py = Axis1[1] * Axis1[3];
251  float NewAxis1_Pz = Axis1[2] * Axis1[3];
252  float NewAxis1_E = Axis1[4];
253 
254  float NewAxis2_Px = Axis2[0] * Axis2[3];
255  float NewAxis2_Py = Axis2[1] * Axis2[3];
256  float NewAxis2_Pz = Axis2[2] * Axis2[3];
257  float NewAxis2_E = Axis2[4];
258 
259  if (Object_Group[i] == 1) {
260  NewAxis1_Px = NewAxis1_Px - (*Object[i]).px();
261  NewAxis1_Py = NewAxis1_Py - (*Object[i]).py();
262  NewAxis1_Pz = NewAxis1_Pz - (*Object[i]).pz();
263  NewAxis1_E = NewAxis1_E - (*Object[i]).energy();
264 
265  } else if (Object_Group[i] == 2) {
266  NewAxis2_Px = NewAxis2_Px - (*Object[i]).px();
267  NewAxis2_Py = NewAxis2_Py - (*Object[i]).py();
268  NewAxis2_Pz = NewAxis2_Pz - (*Object[i]).pz();
269  NewAxis2_E = NewAxis2_E - (*Object[i]).energy();
270  }
271 
272  float mass1 =
273  NewAxis1_E -
274  (((*Object[i]).px() * NewAxis1_Px + (*Object[i]).py() * NewAxis1_Py + (*Object[i]).pz() * NewAxis1_Pz) /
275  (*Object[i]).p());
276 
277  float mass2 =
278  NewAxis2_E -
279  (((*Object[i]).px() * NewAxis2_Px + (*Object[i]).py() * NewAxis2_Py + (*Object[i]).pz() * NewAxis2_Pz) /
280  (*Object[i]).p());
281 
282  if (hemi_meth == 3) {
283  mass1 *= NewAxis1_E / ((NewAxis1_E + (*Object[i]).energy()) * (NewAxis1_E + (*Object[i]).energy()));
284 
285  mass2 *= NewAxis2_E / ((NewAxis2_E + (*Object[i]).energy()) * (NewAxis2_E + (*Object[i]).energy()));
286  }
287 
288  if (mass1 < mass2) {
289  if (Object_Group[i] != 1) {
290  I_Move = true;
291  }
292  Object_Group[i] = 1;
293 
294  Sum1_Px += (*Object[i]).px();
295  Sum1_Py += (*Object[i]).py();
296  Sum1_Pz += (*Object[i]).pz();
297  Sum1_P += (*Object[i]).p();
298  Sum1_E += (*Object[i]).energy();
299  } else {
300  if (Object_Group[i] != 2) {
301  I_Move = true;
302  }
303  Object_Group[i] = 2;
304  Sum2_Px += (*Object[i]).px();
305  Sum2_Py += (*Object[i]).py();
306  Sum2_Pz += (*Object[i]).pz();
307  Sum2_P += (*Object[i]).p();
308  Sum2_E += (*Object[i]).energy();
309  }
310 
311  } else {
312  if (Object_Group[i] == 1) {
313  Sum1_Px += (*Object[i]).px();
314  Sum1_Py += (*Object[i]).py();
315  Sum1_Pz += (*Object[i]).pz();
316  Sum1_P += (*Object[i]).p();
317  Sum1_E += (*Object[i]).energy();
318  }
319  if (Object_Group[i] == 2) {
320  Sum2_Px += (*Object[i]).px();
321  Sum2_Py += (*Object[i]).py();
322  Sum2_Pz += (*Object[i]).pz();
323  Sum2_P += (*Object[i]).p();
324  Sum2_E += (*Object[i]).energy();
325  }
326  }
327  }
328  }
329  } else {
330  throw cms::Exception("Configuration") << "Please give a valid hemisphere association method!";
331  }
332 
333  // recomputing the axes
334 
335  Axis1[3] = sqrt(Sum1_Px * Sum1_Px + Sum1_Py * Sum1_Py + Sum1_Pz * Sum1_Pz);
336  if (Axis1[3] < 0.0001) {
337  edm::LogWarning("HemisphereAlgo") << "ZERO objects in group 1! ";
338  } else {
339  Axis1[0] = Sum1_Px / Axis1[3];
340  Axis1[1] = Sum1_Py / Axis1[3];
341  Axis1[2] = Sum1_Pz / Axis1[3];
342  Axis1[4] = Sum1_E;
343  }
344 
345  Axis2[3] = sqrt(Sum2_Px * Sum2_Px + Sum2_Py * Sum2_Py + Sum2_Pz * Sum2_Pz);
346  if (Axis2[3] < 0.0001) {
347  edm::LogWarning("HemisphereAlgo") << "ZERO objects in group 2!";
348  } else {
349  Axis2[0] = Sum2_Px / Axis2[3];
350  Axis2[1] = Sum2_Py / Axis2[3];
351  Axis2[2] = Sum2_Pz / Axis2[3];
352  Axis2[4] = Sum2_E;
353  }
354 
355  LogDebug("HemisphereAlgo") << " Grouping = ";
356  for (int i = 0; i < vsize; i++) {
357  LogTrace("HemisphereAlgo") << " " << Object_Group[i];
358  }
359  LogTrace("HemisphereAlgo") << std::endl;
360  }
361 
362  status = 1;
363  return 1;
364 }

References Axis1, Axis2, PbPb_ZMuSkimMuonDPG_cff::deltaR, HCALHighEnergyHPDFilter_cfi::energy, PVValHelper::eta, Exception, hemi_meth, mps_fire::i, createfilelist::int, dqmiolumiharvest::j, LogDebug, LogTrace, Object, Object_Group, Object_Noseed, AlCaHLTBitMon_ParallelJobs::p, phi, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, seed_meth, mathSSE::sqrt(), and status.

Referenced by getAxis1(), getAxis2(), and getGrouping().

◆ SetMethod()

void HemisphereAlgo::SetMethod ( int  seed_method,
int  hemisphere_association_method 
)
inline

Definition at line 64 of file HemisphereAlgo.h.

64  {
65  seed_meth = seed_method;
66  hemi_meth = hemisphere_association_method;
67  status = 0;
68  } // sets or overwrites the seed and association methods

References hemi_meth, seed_meth, and status.

◆ SetNoSeed()

void HemisphereAlgo::SetNoSeed ( int  object_number)
inline

Definition at line 70 of file HemisphereAlgo.h.

70  {
71  Object_Noseed[object_number] = 1;
72  status = 0;
73  }

References Object_Noseed, and status.

Member Data Documentation

◆ Axis1

std::vector<float> HemisphereAlgo::Axis1
private

Definition at line 84 of file HemisphereAlgo.h.

Referenced by getAxis1(), and reconstruct().

◆ Axis2

std::vector<float> HemisphereAlgo::Axis2
private

Definition at line 85 of file HemisphereAlgo.h.

Referenced by getAxis2(), and reconstruct().

◆ hemi_meth

int HemisphereAlgo::hemi_meth
private

Definition at line 89 of file HemisphereAlgo.h.

Referenced by reconstruct(), and SetMethod().

◆ Object

std::vector<reco::CandidatePtr> HemisphereAlgo::Object
private

Definition at line 79 of file HemisphereAlgo.h.

Referenced by HemisphereAlgo(), and reconstruct().

◆ Object_Group

std::vector<int> HemisphereAlgo::Object_Group
private

Definition at line 81 of file HemisphereAlgo.h.

Referenced by getGrouping(), and reconstruct().

◆ Object_Noseed

std::vector<int> HemisphereAlgo::Object_Noseed
private

Definition at line 82 of file HemisphereAlgo.h.

Referenced by HemisphereAlgo(), reconstruct(), and SetNoSeed().

◆ seed_meth

int HemisphereAlgo::seed_meth
private

Definition at line 88 of file HemisphereAlgo.h.

Referenced by reconstruct(), and SetMethod().

◆ status

int HemisphereAlgo::status
private

Definition at line 90 of file HemisphereAlgo.h.

Referenced by dirstructure.Comparison::__make_image(), dirstructure.Comparison::__repr__(), getAxis1(), getAxis2(), getGrouping(), reconstruct(), SetMethod(), and SetNoSeed().

mps_fire.i
i
Definition: mps_fire.py:355
HemisphereAlgo::Axis1
std::vector< float > Axis1
Definition: HemisphereAlgo.h:84
multPhiCorr_741_25nsDY_cfi.py
py
Definition: multPhiCorr_741_25nsDY_cfi.py:12
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
HemisphereAlgo::reconstruct
int reconstruct()
Definition: HemisphereAlgo.cc:46
HemisphereAlgo::status
int status
Definition: HemisphereAlgo.h:90
PVValHelper::eta
Definition: PVValidationHelpers.h:69
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
HCALHighEnergyHPDFilter_cfi.energy
energy
Definition: HCALHighEnergyHPDFilter_cfi.py:5
PbPb_ZMuSkimMuonDPG_cff.deltaR
deltaR
Definition: PbPb_ZMuSkimMuonDPG_cff.py:63
edm::LogWarning
Definition: MessageLogger.h:141
HemisphereAlgo::Object_Noseed
std::vector< int > Object_Noseed
Definition: HemisphereAlgo.h:82
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:670
DeltaR
Definition: DeltaR.py:1
HemisphereAlgo::Axis2
std::vector< float > Axis2
Definition: HemisphereAlgo.h:85
createfilelist.int
int
Definition: createfilelist.py:10
HemisphereAlgo::Object_Group
std::vector< int > Object_Group
Definition: HemisphereAlgo.h:81
HemisphereAlgo::hemi_meth
int hemi_meth
Definition: HemisphereAlgo.h:89
DDAxes::phi
multPhiCorr_741_25nsDY_cfi.px
px
Definition: multPhiCorr_741_25nsDY_cfi.py:10
InvariantMass
Definition: InvariantMass.h:6
Exception
Definition: hltDiff.cc:246
HemisphereAlgo::Object
std::vector< reco::CandidatePtr > Object
Definition: HemisphereAlgo.h:79
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:671
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
HemisphereAlgo::seed_meth
int seed_meth
Definition: HemisphereAlgo.h:88
Generated for CMSSW Reference Manual by doxygen 1.8.16