Belle II Software  release-08-01-10
ECLCalDigitVariables.cc
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 
9 /* ECL headers. */
10 #include <ecl/dataobjects/ECLCalDigit.h>
11 #include <ecl/dataobjects/ECLCellIdMapping.h>
12 #include <ecl/dataobjects/ECLDsp.h>
13 #include <ecl/dataobjects/ECLShower.h>
14 #include <ecl/geometry/ECLGeometryPar.h>
15 
16 /* Basf2 headers. */
17 #include <analysis/VariableManager/Manager.h>
18 #include <analysis/dataobjects/Particle.h>
19 #include <framework/datastore/StoreArray.h>
20 #include <framework/datastore/StoreObjPtr.h>
21 #include <framework/geometry/B2Vector3.h>
22 #include <mdst/dataobjects/ECLCluster.h>
23 #include <mdst/dataobjects/MCParticle.h>
24 #include <mdst/dataobjects/Track.h>
25 #include <tracking/dataobjects/ExtHit.h>
26 
27 /* ROOT headers. */
28 #include <Math/VectorUtil.h>
29 
30 using namespace std;
31 
32 namespace Belle2 {
38  namespace ECLCalDigitVariable {
39 
40  // enum with available data types
41  enum varType {
42  energy = 1,
43  time = 2,
44  timeResolution = 3,
45  twoComponentChi2 = 10,
46  twoComponentTotalEnergy = 11,
47  twoComponentHadronEnergy = 12,
48  twoComponentDiodeEnergy = 13,
49  twoComponentSavedChi2_PhotonHadron = 14,
50  twoComponentSavedChi2_PileUpPhoton = 15,
51  twoComponentSavedChi2_PhotonDiode = 16,
52  twoComponentFitType = 17,
53  weight = 18,
54  phi = 20,
55  theta = 21,
56  phiId = 22,
57  thetaId = 23,
58  cellId = 24,
59  mcenergy = 25,
60  usedforenergy = 26,
61  R_geom = 27, // Requires a geometry environment
62  phiOffset = 30,
63  thetaOffset = 31,
64  phiPointing = 32,
65  thetaPointing = 33,
66  twoComponentHadronEnergyFraction = 41,
67  fractionOfShowerEnergy = 42,
68  phiRelativeToShower = 43,
69  thetaRelativeToShower = 44,
70  cosThetaRelativeToShower = 45,
71  rRelativeToShower = 46,
72  };
73 
74  // enum with available center types
75  enum centerType {
76  maxCell = 0,
77  extCell = 1
78  };
79 
80 
82  int getCenterCell(const Particle* particle)
83  {
84  // get maximum cell id for this cluster (using energy and weights)
85  const ECLCluster* cluster = particle->getECLCluster();
86  if (cluster) {
87  int maxCellId = -1;
88  double maxEnergy = -1.;
89 
90  auto clusterDigitRelations = cluster->getRelationsTo<ECLCalDigit>();
91  for (unsigned int ir = 0; ir < clusterDigitRelations.size(); ++ir) {
92  const auto calDigit = clusterDigitRelations.object(ir);
93  const auto weight = clusterDigitRelations.weight(ir);
94 
95  if (calDigit->getEnergy()*weight > maxEnergy) {
96  maxEnergy = calDigit->getEnergy() * weight;
97  maxCellId = calDigit->getCellId();
98  }
99  }
100 
101  return maxCellId;
102  }
103 
104  return -1;
105  }
106 
108  int getExtCell(const Particle* particle)
109  {
110  Const::EDetector myDetID = Const::EDetector::ECL;
111  Const::ChargedStable hypothesis = Const::pion;
112  int pdgCode = abs(hypothesis.getPDGCode());
113 
114  const Track* track = particle->getTrack();
115  if (track) {
116  for (const auto& extHit : track->getRelationsTo<ExtHit>()) {
117  if (abs(extHit.getPdgCode()) != pdgCode) continue;
118  if ((extHit.getDetectorID() != myDetID)) continue;
119  if (extHit.getStatus() != EXT_ENTER) continue;
120 
121  int copyid = extHit.getCopyID();
122 
123  if (copyid == -1) continue;
124  const int cellid = copyid + 1;
125  return cellid;
126  }
127  }
128 
129  return -1;
130  }
131 
133  // This function requires a geometry environment
134  double getCellIdMagnitude(int cellid)
135  {
137  Belle2::B2Vector3D position = geom->GetCrystalPos(cellid - 1);
138  return position.Mag();
139  }
140 
142  std::vector<std::pair<unsigned int, bool>> calculateListOfCrystalEnergyRankAndQuality(ECLShower* shower)
143  {
144  std::vector<std::pair<unsigned int, bool>> listOfCrystalEnergyRankAndQuality;
145  std::vector<std::tuple<double, unsigned int, bool>> energyToSort;
146 
147  RelationVector<ECLCalDigit> relatedDigits = shower->getRelationsTo<ECLCalDigit>();
148 
149  //energyToSort vector is used for sorting digits by calibrated energy
150  for (unsigned int iRel = 0; iRel < relatedDigits.size(); iRel++) {
151 
152  const auto caldigit = relatedDigits.object(iRel);
153  const auto weight = relatedDigits.weight(iRel);
154  bool goodFit = true;
155 
156  //exclude digits without waveforms
157  const double digitChi2 = caldigit->getTwoComponentChi2();
158  if (digitChi2 < 0) goodFit = false;
159 
160  ECLDsp::TwoComponentFitType digitFitType1 = caldigit->getTwoComponentFitType();
161 
162  //exclude digits with poor chi2
163  if (digitFitType1 == ECLDsp::poorChi2) goodFit = false;
164 
165  //exclude digits with diode-crossing fits
166  if (digitFitType1 == ECLDsp::photonDiodeCrossing) goodFit = false;
167 
168  energyToSort.emplace_back(caldigit->getEnergy()*weight, iRel, goodFit);
169  }
170 
171  // sort the vector
172  std::sort(energyToSort.begin(), energyToSort.end(), std::greater<>());
173 
174  for (unsigned int iSorted = 0; iSorted < energyToSort.size(); iSorted++) {
175  listOfCrystalEnergyRankAndQuality.push_back(std::make_pair(std::get<1>(energyToSort[iSorted]),
176  std::get<2>(energyToSort[iSorted])));
177  }
178  return listOfCrystalEnergyRankAndQuality;
179  }
180 
181  ECLShower* getECLShowerFromParticle(const Particle* particle)
182  {
183  const ECLCluster* cluster = particle->getECLCluster();
184  if (!cluster) return nullptr;
185  const auto relShowers = cluster->getRelationsWith<ECLShower>();
186  if (relShowers.size() == 0) return nullptr;
187 
188  if (relShowers.size() == 1) {
189  return relShowers.object(0);
190  } else {
191  B2FATAL("Somehow found more than 1 ECLShower matched to the ECLCluster. This should not be possible!");
192  return nullptr;
193  }
194  }
195 
197  double getCalDigitExpertByEnergyRank(const Particle* particle, const std::vector<double>& vars)
198  {
199 
200  if (!((vars.size() == 2) || (vars.size() == 3))) {
201  B2FATAL("Need two or three parameters (energy rank, variable id, [onlyGoodQualityPSDFits]).");
202  }
203 
204  if (int(std::lround(vars[0])) < 0) B2FATAL("Index cannot be negative.");
205 
206  const unsigned int indexIn = int(std::lround(vars[0]));
207  const int varid = int(std::lround(vars[1]));
208 
209  bool onlyGoodQualityPSDFits = false;
210  if (vars.size() == 3) {
211  onlyGoodQualityPSDFits = static_cast<bool>(std::lround(vars[2]));
212  }
213 
214  ECLShower* shower = getECLShowerFromParticle(particle);
215  if (!shower) return std::numeric_limits<float>::quiet_NaN();
216 
217  // fill the list if it doesn't exist yet.
218  if (shower->getListOfCrystalEnergyRankAndQuality().empty()) {
219  shower->setListOfCrystalEnergyRankAndQuality(calculateListOfCrystalEnergyRankAndQuality(shower));
220  }
221 
222  const std::vector<std::pair<unsigned int, bool>> idxAndQualityList = shower->getListOfCrystalEnergyRankAndQuality();
223 
224  // return nan if we ask for the nth crystal when there are less than n in the shower
225  if (indexIn >= idxAndQualityList.size()) return std::numeric_limits<float>::quiet_NaN();
226 
227  auto relatedDigits = shower->getRelationsTo<ECLCalDigit>();
228  const auto calDigitIndex = idxAndQualityList.at(indexIn).first;
229  const auto goodFit = idxAndQualityList.at(indexIn).second;
230  const auto caldigitSelected = relatedDigits.object(calDigitIndex);
231  const auto weight = relatedDigits.weight(calDigitIndex);
232  const auto digitEnergy = caldigitSelected->getEnergy() * weight;
233 
234  // Mapping object for phi & theta
235  StoreObjPtr<ECLCellIdMapping> mapping;
236  if ((!mapping) and ((varid == varType::phi) or (varid == varType::theta))) {
237  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
238  return std::numeric_limits<double>::quiet_NaN();
239  }
240 
241  // veto bad fits for PSD info if requested
242  if (onlyGoodQualityPSDFits && (!goodFit) && ((varid == varType::twoComponentChi2) ||
243  (varid == varType::twoComponentTotalEnergy) ||
244  (varid == varType::twoComponentHadronEnergy) ||
245  (varid == varType::twoComponentDiodeEnergy) ||
246  (varid == varType::twoComponentFitType) ||
247  (varid == varType::twoComponentHadronEnergyFraction)
248  )) {
249  return std::numeric_limits<double>::quiet_NaN();
250  }
251 
252  if (varid == varType::energy) {
253  return caldigitSelected->getEnergy();
254  } else if (varid == varType::time) {
255  return caldigitSelected->getTime();
256  } else if (varid == varType::twoComponentChi2) {
257  return caldigitSelected->getTwoComponentChi2();
258  } else if (varid == varType::twoComponentTotalEnergy) {
259  return caldigitSelected->getTwoComponentTotalEnergy();
260  } else if (varid == varType::twoComponentHadronEnergy) {
261  return caldigitSelected->getTwoComponentHadronEnergy();
262  } else if (varid == varType::twoComponentSavedChi2_PhotonHadron) {
263  return caldigitSelected->getTwoComponentSavedChi2(ECLDsp::photonHadron);
264  } else if (varid == varType::twoComponentSavedChi2_PileUpPhoton) {
265  return caldigitSelected->getTwoComponentSavedChi2(ECLDsp::photonHadronBackgroundPhoton);
266  } else if (varid == varType::twoComponentSavedChi2_PhotonDiode) {
267  return caldigitSelected->getTwoComponentSavedChi2(ECLDsp::photonDiodeCrossing);
268  } else if (varid == varType::twoComponentDiodeEnergy) {
269  return caldigitSelected->getTwoComponentDiodeEnergy();
270  } else if (varid == varType::twoComponentFitType) {
271  return int(caldigitSelected->getTwoComponentFitType());
272  } else if (varid == varType::cellId) {
273  return caldigitSelected->getCellId();
274  } else if (varid == varType::weight) {
275  return weight;
276  } else if (varid == varType::phi) {
277  return mapping->getCellIdToPhi(caldigitSelected->getCellId());
278  } else if (varid == varType::theta) {
279  return mapping->getCellIdToTheta(caldigitSelected->getCellId());
280  } else if (varid == varType::R_geom) {
281  return getCellIdMagnitude(caldigitSelected->getCellId());
282  } else if (varid == varType::twoComponentHadronEnergyFraction) {
283  if (caldigitSelected-> getTwoComponentTotalEnergy() > 0) {
284  return caldigitSelected->getTwoComponentHadronEnergy() / caldigitSelected->getTwoComponentTotalEnergy();
285  } else {
286  return 0.0;
287  }
288  } else if (varid == varType::fractionOfShowerEnergy) {
289  return digitEnergy / shower->getEnergy();
290 
291  } else if ((varid == varType::phiRelativeToShower) ||
292  (varid == varType::thetaRelativeToShower) ||
293  (varid == varType::cosThetaRelativeToShower) ||
294  (varid == varType::rRelativeToShower)) {
295  ECL::ECLGeometryPar* geometry = ECL::ECLGeometryPar::Instance();
296  B2Vector3D calDigitPosition = geometry->GetCrystalPos(caldigitSelected->getCellId() - 1);
297  B2Vector3D showerPosition;
298  showerPosition.SetMagThetaPhi(shower->getR(), shower->getTheta(), shower->getPhi());
299 
300  ROOT::Math::XYZVector tempP = showerPosition - calDigitPosition;
301  if (varid == varType::rRelativeToShower) return tempP.R();
302  if (varid == varType::thetaRelativeToShower) return tempP.Theta();
303  if (varid == varType::cosThetaRelativeToShower) return tempP.Z() / tempP.R();
304  if (varid == varType::phiRelativeToShower) return tempP.Phi();
305  } else {
306  B2FATAL("variable id not found.");
307  }
308  return std::numeric_limits<double>::quiet_NaN();
309  }
310 
312  double getCalDigitExpert(const Particle* particle, const std::vector<double>& vars)
313  {
314  if (vars.size() != 4) {
315  B2FATAL("Need exactly four parameters (cellid, neighbour area size, variable id, and cluster center (0) or ext (1)).");
316  }
317 
318  StoreObjPtr<ECLCellIdMapping> mapping;
319  const unsigned int posid = int(std::lround(vars[0]));
320 
321  const int nneighbours = int(std::lround(vars[1]));
322  const int varid = int(std::lround(vars[2]));
323  const int extid = int(std::lround(vars[3]));
324 
325  if (!mapping) {
326  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
327  return std::numeric_limits<double>::quiet_NaN();
328  }
329  if (nneighbours != 5 and nneighbours != 7 and nneighbours != 9 and nneighbours != 11) {
330  B2FATAL("Please request 5, 7, 9 or 11 neighbour area.");
331  return std::numeric_limits<double>::quiet_NaN();
332  }
333 
334  // get maximum cell id for this cluster (ignore weights)
335  // TODO: if eclshowers exist we could skip that step.
336  int maxCellId = -1;
337  if (extid == centerType::extCell) {
338  maxCellId = getExtCell(particle);
339  } else {
340  maxCellId = getCenterCell(particle);
341  }
342 
343  if (maxCellId < 0) return std::numeric_limits<double>::quiet_NaN();
344 
345  // get the requested neighbourid
346  int neighbourid = -1;
347  std::vector<short int> neighbours;
348 
349  if (nneighbours == 5) {
350  neighbours = mapping->getCellIdToNeighbour5(maxCellId);
351  } else if (nneighbours == 7) {
352  neighbours = mapping->getCellIdToNeighbour7(maxCellId);
353  } else if (nneighbours == 9) {
354  neighbours = mapping->getCellIdToNeighbour9(maxCellId);
355  } else if (nneighbours == 11) {
356  neighbours = mapping->getCellIdToNeighbour11(maxCellId);
357  }
358 
359  if (posid < neighbours.size()) {
360  neighbourid = neighbours[posid];
361  } else {
362  B2WARNING("This position id is not contained in the requested neighbours.");
363  return std::numeric_limits<double>::quiet_NaN();
364  }
365 
366  // some variables can be returned even if no ECLCalDigit is present
367  if (varid == varType::phi) {
368  return mapping->getCellIdToPhi(neighbourid);
369  } else if (varid == varType::theta) {
370  return mapping->getCellIdToTheta(neighbourid);
371  } else if (varid == varType::R_geom) {
372  return getCellIdMagnitude(neighbourid);
373  } else if (varid == varType::phiId) {
374  return mapping->getCellIdToPhiId(neighbourid);
375  } else if (varid == varType::thetaId) {
376  return mapping->getCellIdToThetaId(neighbourid);
377  } else if (varid == varType::cellId) {
378  return neighbourid;
379  }
380  //... and some really need a ECLCalDigit being present
381  else {
382  const int storearraypos = mapping->getCellIdToStoreArray(neighbourid);
383  StoreArray<ECLCalDigit> eclCalDigits;
384 
385  if (storearraypos >= 0) {
386  if (varid == varType::energy) {
387  return eclCalDigits[storearraypos]->getEnergy();
388  } else if (varid == varType::weight) {
389  const ECLCluster* cluster = particle->getECLCluster();
390  if (cluster == nullptr) {return std::numeric_limits<double>::quiet_NaN();}
391  double weight = 0.;
392  auto relatedDigits = cluster->getRelationsTo<ECLCalDigit>();
393  for (unsigned int iDigit = 0; iDigit < relatedDigits.size(); iDigit++) {
394  const auto caldigit = relatedDigits.object(iDigit);
395  const int digitCellID = caldigit->getCellId();
396  if (digitCellID == neighbourid) {
397  weight = relatedDigits.weight(iDigit);
398  break;
399  }
400  }
401  return weight;
402 
403  } else if (varid == varType::time) {
404  return eclCalDigits[storearraypos]->getTime();
405  } else if (varid == varType::timeResolution) {
406  return eclCalDigits[storearraypos]->getTimeResolution();
407  } else if (varid == varType::twoComponentChi2) {
408  return eclCalDigits[storearraypos]->getTwoComponentChi2();
409  } else if (varid == varType::twoComponentTotalEnergy) {
410  return eclCalDigits[storearraypos]->getTwoComponentTotalEnergy();
411  } else if (varid == varType::twoComponentHadronEnergy) {
412  return eclCalDigits[storearraypos]->getTwoComponentHadronEnergy();
413  } else if (varid == varType::twoComponentSavedChi2_PhotonHadron) {
414  return eclCalDigits[storearraypos]->getTwoComponentSavedChi2(ECLDsp::photonHadron);
415  } else if (varid == varType::twoComponentSavedChi2_PileUpPhoton) {
416  return eclCalDigits[storearraypos]->getTwoComponentSavedChi2(ECLDsp::photonHadronBackgroundPhoton);
417  } else if (varid == varType::twoComponentSavedChi2_PhotonDiode) {
418  return eclCalDigits[storearraypos]->getTwoComponentSavedChi2(ECLDsp::photonDiodeCrossing);
419  } else if (varid == varType::twoComponentDiodeEnergy) {
420  return eclCalDigits[storearraypos]->getTwoComponentDiodeEnergy();
421  } else if (varid == varType::twoComponentFitType) {
422  return int(eclCalDigits[storearraypos]->getTwoComponentFitType());
423  } else if (varid == varType::mcenergy) {
424  // loop over all related MCParticles
425  auto digitMCRelations = eclCalDigits[storearraypos]->getRelationsTo<MCParticle>();
426  double edep = 0.0;
427  for (unsigned int i = 0; i < digitMCRelations.size(); ++i) {
428  edep += digitMCRelations.weight(i);
429  }
430  return edep;
431  } else if (varid == varType::usedforenergy) {
432  const ECLCluster* cluster = particle->getECLCluster();
433  if (cluster) {
434 
435  unsigned int cellid = eclCalDigits[storearraypos]->getCellId();
436  std::vector<unsigned int> listCellIds;
437 
438  auto clusterShowerRelations = cluster->getRelationsWith<ECLShower>();
439 
440  if (clusterShowerRelations.size() == 1) {
441  listCellIds = clusterShowerRelations.object(0)->getListOfCrystalsForEnergy();
442  } else {
443  B2ERROR("Somehow found more than 1 ECLShower matched to the ECLCluster. This should not be possible!");
444  }
445 
446  if (std::find(listCellIds.begin(), listCellIds.end(), cellid) != listCellIds.end()) { //find is faster than count
447  return 1;
448  }
449 
450  return 0;
451  }
452  return std::numeric_limits<double>::quiet_NaN();
453 
454  }
455  } else {
456  return std::numeric_limits<double>::quiet_NaN();
457  }
458  }
459 
460  return std::numeric_limits<double>::quiet_NaN();
461  }
462 
463  double getExtCellExpert(const Particle* particle, int varid, bool front)
464  {
465  ECL::ECLGeometryPar* geometry = ECL::ECLGeometryPar::Instance();
466  if (!geometry) {
467  B2ERROR("Geometry not found!");
468  return std::numeric_limits<double>::quiet_NaN();
469  }
470  const Track* track = particle->getTrack();
471  if (track) {
472  ExtHit* edgeExtHit = nullptr;
473  if (front) {
474  for (const auto& extHit : track->getRelationsTo<ExtHit>()) {
475  if (extHit.getDetectorID() != Const::EDetector::ECL) continue;
476  if (extHit.getStatus() != EXT_ENTER) continue;
477  int crystalID = extHit.getCopyID() - 1;
478  if (crystalID == -1) continue;
479  edgeExtHit = new ExtHit(extHit);
480  break;
481  }
482  } else {
483  auto extHits = track->getRelationsTo<ExtHit>();
484  for (unsigned int iextHit(extHits.size() - 1); iextHit > 0; --iextHit) {
485  const auto extHit = extHits[iextHit];
486  if (extHit->getDetectorID() != Const::EDetector::ECL) continue;
487  if (extHit->getStatus() != EXT_EXIT) continue;
488  int crystalID = extHit->getCopyID() - 1;
489  if (crystalID == -1) break;
490  edgeExtHit = new ExtHit(*extHit);
491  break;
492  }
493  }
494 
495  if (!edgeExtHit) return std::numeric_limits<double>::quiet_NaN();
496  const ROOT::Math::XYZVector& extHitPosition = edgeExtHit->getPosition();
497  const ROOT::Math::XYZVector& trackPointing = edgeExtHit->getMomentum();
498 
499  geometry->Mapping(edgeExtHit->getCopyID() - 1);
500  const int thetaID = geometry->GetThetaID();
501  const int phiID = geometry->GetPhiID();
502  const int cellID = geometry->GetCellID(thetaID, phiID);
503 
504  const ROOT::Math::XYZVector& crystalCenterPosition =
505  geometry->GetCrystalPos(cellID);
506  const ROOT::Math::XYZVector& crystalOrientation =
507  geometry->GetCrystalVec(cellID);
508  const ROOT::Math::XYZVector& crystalPositionOnSurface =
509  crystalCenterPosition -
510  (crystalCenterPosition - extHitPosition).Dot(
511  crystalOrientation.Unit()) * crystalOrientation.Unit();
512  if (varid == varType::phiOffset) {
513  return ROOT::Math::VectorUtil::Phi_mpi_pi(extHitPosition.Phi() - crystalPositionOnSurface.Phi());
514  } else if (varid == varType::thetaOffset) {
515  return extHitPosition.Theta() - crystalPositionOnSurface.Theta();
516  } else if (varid == varType::phiPointing) {
517  return ROOT::Math::VectorUtil::Phi_mpi_pi(trackPointing.Phi() - crystalOrientation.Phi());
518  } else if (varid == varType::thetaPointing) {
519  return trackPointing.Theta() - crystalOrientation.Theta();
520  }
521  }
522 
523  return std::numeric_limits<double>::quiet_NaN();
524  }
525  }
526 
527  namespace Variable {
528 
530  double getECLCalDigitEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
531  {
532 
533  if (vars.size() != 1) {
534  B2FATAL("Need exactly one parameter (energy index).");
535  }
536  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::energy};
537  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
538  }
539 
541  double getECLCalDigitTimeByEnergyRank(const Particle* particle, const std::vector<double>& vars)
542  {
543  if (vars.size() != 1) {
544  B2FATAL("Need exactly one parameters (energy index).");
545  }
546  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::time};
547  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
548  }
549 
551  double getECLCalDigitMCEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
552  {
553  if (vars.size() != 1) {
554  B2FATAL("Need exactly one parameters (energy index).");
555  }
556  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::mcenergy};
557  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
558  }
559 
561  double getCellIdByEnergyRank(const Particle* particle, const std::vector<double>& vars)
562  {
563  if (vars.size() != 1) {
564  B2FATAL("Need exactly one parameters (energy index).");
565  }
566  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::cellId};
567  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
568  }
569 
571  double getTwoComponentFitTypeByEnergyRank(const Particle* particle, const std::vector<double>& vars)
572  {
573  if (!((vars.size() == 1) || (vars.size() == 2))) {
574  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
575  }
576  double onlyGoodQualityPSDFits = 0.0;
577  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
578 
579  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentFitType, onlyGoodQualityPSDFits};
580  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
581  }
582 
584  double getTwoComponentChi2ByEnergyRank(const Particle* particle, const std::vector<double>& vars)
585  {
586  if (!((vars.size() == 1) || (vars.size() == 2))) {
587  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
588  }
589  double onlyGoodQualityPSDFits = 0.0;
590  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
591 
592  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentChi2, onlyGoodQualityPSDFits};
593  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
594  }
595 
597  double getTwoComponentTotalEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
598  {
599  if (!((vars.size() == 1) || (vars.size() == 2))) {
600  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
601  }
602  double onlyGoodQualityPSDFits = 0.0;
603  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
604 
605  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentTotalEnergy, onlyGoodQualityPSDFits};
606  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
607  }
608 
610  double getTwoComponentHadronEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
611  {
612  if (!((vars.size() == 1) || (vars.size() == 2))) {
613  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
614  }
615  double onlyGoodQualityPSDFits = 0.0;
616  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
617 
618  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentHadronEnergy, onlyGoodQualityPSDFits};
619  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
620  }
621 
623  double getTwoComponentHadronEnergyFractionByEnergyRank(const Particle* particle, const std::vector<double>& vars)
624  {
625  if (!((vars.size() == 1) || (vars.size() == 2))) {
626  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
627  }
628  double onlyGoodQualityPSDFits = 0.0;
629  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
630 
631  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentHadronEnergyFraction, onlyGoodQualityPSDFits};
632  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
633  }
634 
636  double getTwoComponentDiodeEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
637  {
638  if (!((vars.size() == 1) || (vars.size() == 2))) {
639  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
640  }
641  double onlyGoodQualityPSDFits = 0.0;
642  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
643 
644  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentDiodeEnergy, onlyGoodQualityPSDFits};
645  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
646  }
647 
649  double getTwoComponentChi2SavedByEnergyRank_PhotonHadron(const Particle* particle, const std::vector<double>& vars)
650  {
651  if (!((vars.size() == 1) || (vars.size() == 2))) {
652  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
653  }
654  double onlyGoodQualityPSDFits = 0.0;
655  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
656 
657  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonHadron, onlyGoodQualityPSDFits};
658  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
659  }
660 
662  double getTwoComponentChi2SavedByEnergyRank_PileUpPhoton(const Particle* particle, const std::vector<double>& vars)
663  {
664  if (!((vars.size() == 1) || (vars.size() == 2))) {
665  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
666  }
667  double onlyGoodQualityPSDFits = 0.0;
668  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
669 
670  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentSavedChi2_PileUpPhoton, onlyGoodQualityPSDFits};
671  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
672  }
673 
675  double getTwoComponentChi2SavedByEnergyRank_PhotonDiode(const Particle* particle, const std::vector<double>& vars)
676  {
677  if (!((vars.size() == 1) || (vars.size() == 2))) {
678  B2FATAL("Need one or two parameters (energy index, [onlyGoodQualityPSDFits]).");
679  }
680  double onlyGoodQualityPSDFits = 0.0;
681  if (vars.size() == 2) {onlyGoodQualityPSDFits = std::lround(vars[1]);}
682 
683  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonDiode, onlyGoodQualityPSDFits};
684  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
685  }
686 
688  double getWeightByEnergyRank(const Particle* particle, const std::vector<double>& vars)
689  {
690  if (vars.size() != 1) {
691  B2FATAL("Need exactly one parameters (energy index).");
692  }
693  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::weight};
694  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
695  }
696 
698  double getECLCalDigitFractionOfShowerEnergyByEnergyRank(const Particle* particle, const std::vector<double>& vars)
699  {
700  if (vars.size() != 1) {
701  B2FATAL("Need exactly one parameters (energy index).");
702  }
703  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::fractionOfShowerEnergy};
704  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
705  }
706 
708  double getECLCalDigitPhiRelativeToShowerByEnergyRank(const Particle* particle, const std::vector<double>& vars)
709  {
710  if (vars.size() != 1) {
711  B2FATAL("Need exactly one parameters (energy index).");
712  }
713  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::phiRelativeToShower};
714  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
715  }
716 
718  double getECLCalDigitThetaRelativeToShowerByEnergyRank(const Particle* particle, const std::vector<double>& vars)
719  {
720  if (vars.size() != 1) {
721  B2FATAL("Need exactly one parameters (energy index).");
722  }
723  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::thetaRelativeToShower};
724  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
725  }
726 
728  double getECLCalDigitCosThetaRelativeToShowerByEnergyRank(const Particle* particle, const std::vector<double>& vars)
729  {
730  if (vars.size() != 1) {
731  B2FATAL("Need exactly one parameters (energy index).");
732  }
733  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::cosThetaRelativeToShower};
734  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
735  }
736 
738  double getECLCalDigitRadiusRelativeToShowerByEnergyRank(const Particle* particle, const std::vector<double>& vars)
739  {
740  if (vars.size() != 1) {
741  B2FATAL("Need exactly one parameters (energy index).");
742  }
743  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::rRelativeToShower};
744  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
745  }
746 
748  double getECLCalDigitEnergy(const Particle* particle, const std::vector<double>& vars)
749  {
750  if (vars.size() != 2) {
751  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
752  }
753 
754  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::energy, ECLCalDigitVariable::centerType::maxCell};
755  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
756  }
757 
758 
760  double getECLCalDigitWeight(const Particle* particle, const std::vector<double>& vars)
761  {
762  if (vars.size() != 2) {
763  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
764  }
765 
766  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::weight, ECLCalDigitVariable::centerType::maxCell};
767  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
768  }
769 
771  double getMCEnergy(const Particle* particle, const std::vector<double>& vars)
772  {
773  if (vars.size() != 2) {
774  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
775  }
776 
777  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::mcenergy, ECLCalDigitVariable::centerType::maxCell};
778  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
779  }
780 
781 
783  double getExtECLCalDigitEnergy(const Particle* particle, const std::vector<double>& vars)
784  {
785  if (vars.size() != 2) {
786  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
787  }
788 
789  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::energy, ECLCalDigitVariable::centerType::extCell};
790  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
791  }
792 
793 
795  double getECLCalDigitTime(const Particle* particle, const std::vector<double>& vars)
796  {
797  if (vars.size() != 2) {
798  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
799  }
800 
801  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::time, ECLCalDigitVariable::centerType::maxCell};
802  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
803  }
804 
806  double getExtECLCalDigitTime(const Particle* particle, const std::vector<double>& vars)
807  {
808  if (vars.size() != 2) {
809  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
810  }
811 
812  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::time, ECLCalDigitVariable::centerType::extCell};
813  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
814  }
815 
817  double getExtECLCalDigitTimeResolution(const Particle* particle, const std::vector<double>& vars)
818  {
819  if (vars.size() != 2) {
820  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
821  }
822 
823  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::timeResolution, ECLCalDigitVariable::centerType::extCell};
824  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
825  }
826 
828  double getECLCalDigitTimeResolution(const Particle* particle, const std::vector<double>& vars)
829  {
830  if (vars.size() != 2) {
831  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
832  }
833 
834  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::timeResolution, ECLCalDigitVariable::centerType::maxCell};
835  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
836  }
837 
839  double getTwoComponentChi2(const Particle* particle, const std::vector<double>& vars)
840  {
841  if (vars.size() != 2) {
842  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
843  }
844 
845  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentChi2, ECLCalDigitVariable::centerType::maxCell};
846  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
847  }
848 
850  double getExtECLCalDigitTwoComponentChi2(const Particle* particle, const std::vector<double>& vars)
851  {
852  if (vars.size() != 2) {
853  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
854  }
855 
856  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentChi2, ECLCalDigitVariable::centerType::extCell};
857  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
858  }
859 
861  double getTwoComponentTotalEnergy(const Particle* particle, const std::vector<double>& vars)
862  {
863  if (vars.size() != 2) {
864  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
865  }
866 
867  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentTotalEnergy, ECLCalDigitVariable::centerType::maxCell};
868  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
869  }
870 
872  double getExtECLCalDigitTwoComponentTotalEnergy(const Particle* particle, const std::vector<double>& vars)
873  {
874  if (vars.size() != 2) {
875  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
876  }
877 
878  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentTotalEnergy, ECLCalDigitVariable::centerType::extCell};
879  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
880  }
881 
883  double getTwoComponentHadronEnergy(const Particle* particle, const std::vector<double>& vars)
884  {
885  if (vars.size() != 2) {
886  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
887  }
888 
889  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentHadronEnergy, ECLCalDigitVariable::centerType::maxCell};
890  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
891  }
892 
894  double getUsedForClusterEnergy(const Particle* particle, const std::vector<double>& vars)
895  {
896  if (vars.size() != 2) {
897  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
898  }
899 
900  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::usedforenergy, ECLCalDigitVariable::centerType::maxCell};
901  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
902  }
903 
905  double getExtECLCalDigitTwoComponentHadronEnergy(const Particle* particle, const std::vector<double>& vars)
906  {
907  if (vars.size() != 2) {
908  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
909  }
910 
911  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentHadronEnergy, ECLCalDigitVariable::centerType::extCell};
912  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
913  }
914 
916  double getECLCalDigitTwoComponentDiodeEnergy(const Particle* particle, const std::vector<double>& vars)
917  {
918  if (vars.size() != 2) {
919  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
920  }
921 
922  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentDiodeEnergy, ECLCalDigitVariable::centerType::maxCell};
923  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
924  }
925 
927  double getExtECLCalDigitTwoComponentDiodeEnergy(const Particle* particle, const std::vector<double>& vars)
928  {
929  if (vars.size() != 2) {
930  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
931  }
932 
933  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentDiodeEnergy, ECLCalDigitVariable::centerType::extCell};
934  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
935  }
936 
938  double getECLCalDigitTwoComponentFitType(const Particle* particle, const std::vector<double>& vars)
939  {
940  if (vars.size() != 2) {
941  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
942  }
943 
944  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentFitType, ECLCalDigitVariable::centerType::maxCell};
945  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
946  }
947 
949  double getExtECLCalDigitTwoComponentFitType(const Particle* particle, const std::vector<double>& vars)
950  {
951  if (vars.size() != 2) {
952  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
953  }
954 
955  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentFitType, ECLCalDigitVariable::centerType::extCell};
956  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
957  }
958 
960  double getECLCalDigitTwoComponentChi2Saved_PhotonHadron(const Particle* particle, const std::vector<double>& vars)
961  {
962  if (vars.size() != 2) {
963  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
964  }
965 
966  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonHadron, ECLCalDigitVariable::centerType::maxCell};
967  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
968  }
969 
971  double getExtECLCalDigitTwoComponentChi2Saved_PhotonHadron(const Particle* particle, const std::vector<double>& vars)
972  {
973  if (vars.size() != 2) {
974  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
975  }
976 
977  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonHadron, ECLCalDigitVariable::centerType::extCell};
978  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
979  }
980 
982  double getECLCalDigitTwoComponentChi2Saved_PileUpPhoton(const Particle* particle, const std::vector<double>& vars)
983  {
984  if (vars.size() != 2) {
985  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
986  }
987 
988  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PileUpPhoton, ECLCalDigitVariable::centerType::maxCell};
989  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
990  }
991 
993  double getExtECLCalDigitTwoComponentChi2Saved_PileUpPhoton(const Particle* particle, const std::vector<double>& vars)
994  {
995  if (vars.size() != 2) {
996  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
997  }
998 
999  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PileUpPhoton, ECLCalDigitVariable::centerType::extCell};
1000  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1001  }
1002 
1004  double getECLCalDigitTwoComponentChi2Saved_PhotonDiode(const Particle* particle, const std::vector<double>& vars)
1005  {
1006  if (vars.size() != 2) {
1007  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1008  }
1009 
1010  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonDiode, ECLCalDigitVariable::centerType::maxCell};
1011  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1012  }
1013 
1015  double getExtECLCalDigitTwoComponentChi2Saved_PhotonDiode(const Particle* particle, const std::vector<double>& vars)
1016  {
1017  if (vars.size() != 2) {
1018  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1019  }
1020 
1021  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::twoComponentSavedChi2_PhotonDiode, ECLCalDigitVariable::centerType::extCell};
1022  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1023  }
1024 
1026  double getPhi(const Particle* particle, const std::vector<double>& vars)
1027  {
1028  if (vars.size() != 2) {
1029  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1030  }
1031 
1032  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::phi, ECLCalDigitVariable::centerType::maxCell};
1033  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1034  }
1035 
1037  double getExtPhi(const Particle* particle, const std::vector<double>& vars)
1038  {
1039  if (vars.size() != 2) {
1040  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1041  }
1042 
1043  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::phi, ECLCalDigitVariable::centerType::extCell};
1044  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1045  }
1046 
1048  double getTheta(const Particle* particle, const std::vector<double>& vars)
1049  {
1050  if (vars.size() != 2) {
1051  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1052  }
1053 
1054  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::theta, ECLCalDigitVariable::centerType::maxCell};
1055  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1056  }
1057 
1059  double getExtTheta(const Particle* particle, const std::vector<double>& vars)
1060  {
1061  if (vars.size() != 2) {
1062  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1063  }
1064 
1065  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::theta, ECLCalDigitVariable::centerType::extCell};
1066  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1067  }
1068 
1070  double getPhiId(const Particle* particle, const std::vector<double>& vars)
1071  {
1072  if (vars.size() != 2) {
1073  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1074  }
1075 
1076  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::phiId, ECLCalDigitVariable::centerType::maxCell};
1077  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1078  }
1079 
1081  double getExtPhiId(const Particle* particle, const std::vector<double>& vars)
1082  {
1083  if (vars.size() != 2) {
1084  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1085  }
1086 
1087  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::phiId, ECLCalDigitVariable::centerType::extCell};
1088  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1089  }
1090 
1092  double getThetaId(const Particle* particle, const std::vector<double>& vars)
1093  {
1094  if (vars.size() != 2) {
1095  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1096  }
1097 
1098  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::thetaId, ECLCalDigitVariable::centerType::maxCell};
1099  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1100  }
1101 
1103  double getExtThetaId(const Particle* particle, const std::vector<double>& vars)
1104  {
1105  if (vars.size() != 2) {
1106  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1107  }
1108 
1109  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::thetaId, ECLCalDigitVariable::centerType::extCell};
1110  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1111  }
1112 
1114  double getCellId(const Particle* particle, const std::vector<double>& vars)
1115  {
1116  if (vars.size() != 2) {
1117  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1118  }
1119 
1120  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::cellId, ECLCalDigitVariable::centerType::maxCell};
1121  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1122  }
1123 
1125  double getCenterCellId(const Particle* particle)
1126  {
1127  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1128 
1129  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1130  return centercellid;
1131  }
1132 
1134  double getCenterCellCrystalTheta(const Particle* particle)
1135  {
1136  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1137  StoreObjPtr<ECLCellIdMapping> mapping;
1138 
1139  if (!mapping) {
1140  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1141  return std::numeric_limits<double>::quiet_NaN();
1142  }
1143 
1144  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1145  return mapping->getCellIdToTheta(centercellid);
1146  }
1147 
1149  double getCenterCellCrystalPhi(const Particle* particle)
1150  {
1151  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1152  StoreObjPtr<ECLCellIdMapping> mapping;
1153 
1154  if (!mapping) {
1155  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1156  return std::numeric_limits<double>::quiet_NaN();
1157  }
1158 
1159  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1160  return mapping->getCellIdToPhi(centercellid);
1161  }
1162 
1164  double getExtCellId(const Particle* particle)
1165  {
1166  const int extcellid = ECLCalDigitVariable::getExtCell(particle);
1167 
1168  if (extcellid < 0) return std::numeric_limits<double>::quiet_NaN();
1169  return extcellid;
1170  }
1171 
1173  double getCenterCellThetaId(const Particle* particle)
1174  {
1175  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1176  StoreObjPtr<ECLCellIdMapping> mapping;
1177 
1178  if (!mapping) {
1179  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1180  return std::numeric_limits<double>::quiet_NaN();
1181  }
1182 
1183  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1184  return mapping->getCellIdToThetaId(centercellid);
1185  }
1186 
1188  double getCenterCellPhiId(const Particle* particle)
1189  {
1190  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1191  StoreObjPtr<ECLCellIdMapping> mapping;
1192 
1193  if (!mapping) {
1194  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1195  return std::numeric_limits<double>::quiet_NaN();
1196  }
1197 
1198  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1199  return mapping->getCellIdToPhiId(centercellid);
1200  }
1201 
1203  double getExtCellThetaId(const Particle* particle)
1204  {
1205  const int extcellid = ECLCalDigitVariable::getExtCell(particle);
1206  StoreObjPtr<ECLCellIdMapping> mapping;
1207 
1208  if (!mapping) {
1209  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1210  return std::numeric_limits<double>::quiet_NaN();
1211  }
1212 
1213  if (extcellid < 0) return std::numeric_limits<double>::quiet_NaN();
1214  return mapping->getCellIdToThetaId(extcellid);
1215  }
1216 
1218  double getExtCellPhiId(const Particle* particle)
1219  {
1220  const int extcellid = ECLCalDigitVariable::getExtCell(particle);
1221  StoreObjPtr<ECLCellIdMapping> mapping;
1222 
1223  if (!mapping) {
1224  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1225  return std::numeric_limits<double>::quiet_NaN();
1226  }
1227 
1228  if (extcellid < 0) return std::numeric_limits<double>::quiet_NaN();
1229  return mapping->getCellIdToPhiId(extcellid);
1230  }
1231 
1233  double getExtCellCrystalTheta(const Particle* particle)
1234  {
1235  const int extcellid = ECLCalDigitVariable::getExtCell(particle);
1236  StoreObjPtr<ECLCellIdMapping> mapping;
1237 
1238  if (!mapping) {
1239  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1240  return std::numeric_limits<double>::quiet_NaN();
1241  }
1242 
1243  if (extcellid < 0) return std::numeric_limits<double>::quiet_NaN();
1244  return mapping->getCellIdToTheta(extcellid);
1245  }
1246 
1248  double getExtCellCrystalPhi(const Particle* particle)
1249  {
1250  const int extcellid = ECLCalDigitVariable::getExtCell(particle);
1251  StoreObjPtr<ECLCellIdMapping> mapping;
1252 
1253  if (!mapping) {
1254  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1255  return std::numeric_limits<double>::quiet_NaN();
1256  }
1257 
1258  if (extcellid < 0) return std::numeric_limits<double>::quiet_NaN();
1259  return mapping->getCellIdToPhi(extcellid);
1260  }
1261 
1263  double getCenterCellIndex(const Particle* particle, const std::vector<double>& vars)
1264  {
1265  if (vars.size() != 1) {
1266  B2FATAL("Need exactly one parameters (neighbour area size).");
1267  }
1268 
1269  StoreObjPtr<ECLCellIdMapping> mapping;
1270  const int nneighbours = int(std::lround(vars[0]));
1271 
1272  if (!mapping) {
1273  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1274  return std::numeric_limits<double>::quiet_NaN();
1275  }
1276  if (nneighbours != 5 and nneighbours != 7 and nneighbours != 9 and nneighbours != 11) {
1277  B2FATAL("Please request 5, 7, 9 or 11 neighbour area.");
1278  return std::numeric_limits<double>::quiet_NaN();
1279  }
1280 
1281  const int centercellid = ECLCalDigitVariable::getCenterCell(particle);
1282  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1283 
1284  std::vector<short int> neighbours;
1285 
1286  if (nneighbours == 5) {
1287  neighbours = mapping->getCellIdToNeighbour5(centercellid);
1288  } else if (nneighbours == 7) {
1289  neighbours = mapping->getCellIdToNeighbour7(centercellid);
1290  } else if (nneighbours == 9) {
1291  neighbours = mapping->getCellIdToNeighbour9(centercellid);
1292  } else if (nneighbours == 11) {
1293  neighbours = mapping->getCellIdToNeighbour11(centercellid);
1294  }
1295 
1296  for (unsigned int idx = 0; idx < neighbours.size(); idx++) {
1297  if (neighbours[idx] == centercellid) return idx;
1298  }
1299 
1300  return std::numeric_limits<double>::quiet_NaN();
1301  }
1302 
1303 
1305  double getExtCenterCellIndex(const Particle* particle, const std::vector<double>& vars)
1306  {
1307  if (vars.size() != 1) {
1308  B2FATAL("Need exactly one parameters (neighbour area size).");
1309  }
1310 
1311  StoreObjPtr<ECLCellIdMapping> mapping;
1312  const int nneighbours = int(std::lround(vars[0]));
1313 
1314  if (!mapping) {
1315  B2ERROR("Mapping not found, did you forget to run the eclFillCellIdMapping module?");
1316  return std::numeric_limits<double>::quiet_NaN();
1317  }
1318  if (nneighbours != 5 and nneighbours != 7 and nneighbours != 9 and nneighbours != 11) {
1319  B2FATAL("Please request 5, 7, 9 or 11 neighbour area.");
1320  return std::numeric_limits<double>::quiet_NaN();
1321  }
1322 
1323  const int centercellid = ECLCalDigitVariable::getExtCell(particle);
1324  if (centercellid < 0) return std::numeric_limits<double>::quiet_NaN();
1325 
1326  std::vector<short int> neighbours;
1327 
1328  if (nneighbours == 5) {
1329  neighbours = mapping->getCellIdToNeighbour5(centercellid);
1330  } else if (nneighbours == 7) {
1331  neighbours = mapping->getCellIdToNeighbour7(centercellid);
1332  } else if (nneighbours == 9) {
1333  neighbours = mapping->getCellIdToNeighbour9(centercellid);
1334  } else if (nneighbours == 11) {
1335  neighbours = mapping->getCellIdToNeighbour11(centercellid);
1336  }
1337 
1338  for (unsigned int idx = 0; idx < neighbours.size(); idx++) {
1339  if (neighbours[idx] == centercellid) return idx;
1340  }
1341 
1342  return std::numeric_limits<double>::quiet_NaN();
1343  }
1344 
1345  double getTotalECLCalDigitMCEnergy(const Particle* particle)
1346  {
1347  const MCParticle* mcparticle = particle->getRelatedTo<MCParticle>();
1348  if (mcparticle == nullptr)
1349  return std::numeric_limits<double>::quiet_NaN();
1350 
1351  double sum = 0.0;
1352  auto mcDigitRelations = mcparticle->getRelationsWith<ECLCalDigit>();
1353  for (unsigned int ir = 0; ir < mcDigitRelations.size(); ++ir) {
1354  sum += mcDigitRelations.weight(ir);
1355  }
1356 
1357  return sum;
1358  }
1359 
1360  double getClusterTotalECLCalDigitMCEnergy(const Particle* particle)
1361  {
1362  const ECLCluster* cluster = particle->getECLCluster();
1363  if (cluster == nullptr) {return std::numeric_limits<double>::quiet_NaN();}
1364  const MCParticle* mcparticle = particle->getRelatedTo<MCParticle>();
1365  if (mcparticle == nullptr) {return std::numeric_limits<double>::quiet_NaN();}
1366 
1367  double sum = 0.0;
1368  auto relatedDigits = cluster->getRelationsTo<ECLCalDigit>();
1369  for (unsigned int iDigit = 0; iDigit < relatedDigits.size(); iDigit++) {
1370  const auto caldigit = relatedDigits.object(iDigit);
1371  auto mcDigitRelations = caldigit->getRelationsTo<MCParticle>();
1372  for (unsigned int i = 0; i < mcDigitRelations.size(); i++) {
1373  const MCParticle* digitmcparticle = mcDigitRelations.object(i);
1374  if (digitmcparticle == mcparticle) {
1375  sum += mcDigitRelations.weight(i);
1376  }
1377  }
1378  }
1379  return sum;
1380 
1381  }
1382 
1383  double getClusterECLCalDigitMCEnergy(const Particle* particle)
1384  {
1385  // get MCParticle (return if there is none)
1386  const MCParticle* mcparticle = particle->getRelatedTo<MCParticle>();
1387  if (mcparticle == nullptr)
1388  return std::numeric_limits<double>::quiet_NaN();
1389 
1390  const ECLCluster* cluster = particle->getECLCluster();
1391  if (cluster) {
1392  std::vector<unsigned int> listCellIds;
1393  auto clusterShowerRelations = cluster->getRelationsWith<ECLShower>();
1394  if (clusterShowerRelations.size() == 1) {
1395  listCellIds = clusterShowerRelations.object(0)->getListOfCrystalsForEnergy();
1396  } else {
1397  B2ERROR("Somehow found more than 1 ECLShower matched to the ECLCluster. This should not be possible!");
1398  }
1399 
1400  double sum = 0.0;
1401  auto clusterDigitRelations = mcparticle->getRelationsWith<ECLCalDigit>();
1402  for (unsigned int ir = 0; ir < clusterDigitRelations.size(); ++ir) {
1403 
1404  // check if this digit is in the current cluster
1405  unsigned int cellid = clusterDigitRelations.object(ir)->getCellId();
1406  if (std::find(listCellIds.begin(), listCellIds.end(), cellid) != listCellIds.end()) { //find is faster than count
1407  sum += clusterDigitRelations.weight(ir);
1408  }
1409  }
1410 
1411  return sum;
1412  }
1413  return std::numeric_limits<float>::quiet_NaN();
1414  }
1415 
1416  double getExtFrontPositionPhiOffset(const Particle* particle)
1417  {
1418  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::phiOffset, true);
1419  }
1420 
1421  double getExtFrontPositionThetaOffset(const Particle* particle)
1422  {
1423  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::thetaOffset, true);
1424  }
1425 
1426  double getExtFrontPositionPhiPointing(const Particle* particle)
1427  {
1428  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::phiPointing, true);
1429  }
1430 
1431  double getExtFrontPositionThetaPointing(const Particle* particle)
1432  {
1433  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::thetaPointing, true);
1434  }
1435 
1436  double getExtBackPositionPhiOffset(const Particle* particle)
1437  {
1438  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::phiOffset, false);
1439  }
1440 
1441  double getExtBackPositionThetaOffset(const Particle* particle)
1442  {
1443  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::thetaOffset, false);
1444  }
1445 
1446  double getExtBackPositionPhiPointing(const Particle* particle)
1447  {
1448  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::phiPointing, false);
1449  }
1450 
1451  double getExtBackPositionThetaPointing(const Particle* particle)
1452  {
1453  return ECLCalDigitVariable::getExtCellExpert(particle, ECLCalDigitVariable::varType::thetaPointing, false);
1454  }
1455 
1457  double getPhiByEnergyRank(const Particle* particle, const std::vector<double>& vars)
1458  {
1459  if (vars.size() != 1) {
1460  B2FATAL("Need exactly one parameters (energy index).");
1461  }
1462  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::phi};
1463  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
1464  }
1465 
1467  double getThetaByEnergyRank(const Particle* particle, const std::vector<double>& vars)
1468  {
1469  if (vars.size() != 1) {
1470  B2FATAL("Need exactly one parameters (energy index).");
1471  }
1472  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::theta};
1473  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
1474  }
1475 
1477  double getR(const Particle* particle, const std::vector<double>& vars)
1478  {
1479  if (vars.size() != 2) {
1480  B2FATAL("Need exactly two parameters (cellid and neighbour area size).");
1481  }
1482 
1483  std::vector<double> parameters {vars[0], vars[1], ECLCalDigitVariable::varType::R_geom, ECLCalDigitVariable::centerType::maxCell};
1484  return ECLCalDigitVariable::getCalDigitExpert(particle, parameters);
1485  }
1486 
1488  double getRByEnergyRank(const Particle* particle, const std::vector<double>& vars)
1489  {
1490  if (vars.size() != 1) {
1491  B2FATAL("Need exactly one parameters (energy index).");
1492  }
1493  std::vector<double> parameters {vars[0], ECLCalDigitVariable::varType::R_geom};
1494  return ECLCalDigitVariable::getCalDigitExpertByEnergyRank(particle, parameters);
1495  }
1496 
1497  double getClusterNHitsThreshold(const Particle* particle, const std::vector<double>& vars)
1498  {
1499  if (vars.size() != 1) {
1500  B2FATAL("Need exactly one parameter (energy threshold in GeV).");
1501  }
1502  const double threshold = vars[0];
1503 
1504  const ECLCluster* cluster = particle->getECLCluster();
1505  if (cluster) {
1506  double nhits = 0.;
1507 
1508  auto clusterDigitRelations = cluster->getRelationsTo<ECLCalDigit>();
1509  for (unsigned int ir = 0; ir < clusterDigitRelations.size(); ++ir) {
1510  const auto calDigit = clusterDigitRelations.object(ir);
1511  const auto weight = clusterDigitRelations.weight(ir);
1512 
1513  // take the unweighted eclcaldigit energy for this check (closer to real hardware threshold)
1514  if (calDigit->getEnergy() > threshold) {
1515  nhits += weight;
1516  }
1517  }
1518 
1519  return nhits;
1520  }
1521  return std::numeric_limits<float>::quiet_NaN();
1522  }
1523 
1524 
1525  VARIABLE_GROUP("ECL Calibration (cDST)");
1526  REGISTER_VARIABLE("eclcaldigitEnergy(i, j)", getECLCalDigitEnergy,
1527  "[calibration] Returns the energy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1528  REGISTER_VARIABLE("eclcaldigitWeight(i, j)", getECLCalDigitWeight,
1529  "[calibration] Returns the weight of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1530  REGISTER_VARIABLE("eclcaldigitTime(i, j)", getECLCalDigitTime,
1531  "[calibration] Returns the time of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1532  REGISTER_VARIABLE("eclcaldigitTimeResolution(i, j)", getECLCalDigitTimeResolution,
1533  "[calibration] Returns the time resolution (dt99) of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1534  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2(i, j)", getTwoComponentChi2,
1535  "[calibration] Returns the two component fit chi2 of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1536  REGISTER_VARIABLE("eclcaldigitTwoComponentTotalEnergy(i, j)", getTwoComponentTotalEnergy,
1537  "[calibration] Returns the two component total energy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1538  REGISTER_VARIABLE("eclcaldigitTwoComponentHadronEnergy(i, j)", getTwoComponentHadronEnergy,
1539  "[calibration] Returns the two component hadron energy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1540  REGISTER_VARIABLE("eclcaldigitPhi(i, j)", getPhi,
1541  "[calibration] Returns phi of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1542  REGISTER_VARIABLE("eclcaldigitTheta(i, j)", getTheta,
1543  "[calibration] Returns theta of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1544  REGISTER_VARIABLE("eclcaldigitR(i, j)", getR,
1545  "Returns R (from a geometry object) of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1546  REGISTER_VARIABLE("eclcaldigitPhiId(i, j)", getPhiId,
1547  "[calibration] Returns the phi Id of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1548  REGISTER_VARIABLE("eclcaldigitThetaId(i, j)", getThetaId,
1549  "[calibration] Returns the theta Id of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1550  REGISTER_VARIABLE("eclcaldigitCellId(i, j)", getCellId,
1551  "[calibration] Returns the cell id of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours (1-based)");
1552  REGISTER_VARIABLE("eclcaldigitUsedForClusterEnergy(i, j)", getUsedForClusterEnergy,
1553  " [calibration] Returns the 0 (not used) 1 (used) of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours (1-based)");
1554 
1555  REGISTER_VARIABLE("eclcaldigitCenterCellId", getCenterCellId, "[calibration] Returns the center cell id");
1556  REGISTER_VARIABLE("eclcaldigitCenterCellThetaId", getCenterCellThetaId, "[calibration] Returns the center cell theta id");
1557  REGISTER_VARIABLE("eclcaldigitCenterCellPhiId", getCenterCellPhiId, "[calibration] Returns the center cell phi id");
1558  REGISTER_VARIABLE("eclcaldigitCenterCellCrystalTheta", getCenterCellCrystalTheta,
1559  "[calibration] Returns the center cell crystal theta");
1560  REGISTER_VARIABLE("eclcaldigitCenterCellCrystalPhi", getCenterCellCrystalPhi,
1561  "[calibration] Returns the center cell crystal phi");
1562  REGISTER_VARIABLE("eclcaldigitCenterCellIndex(i)", getCenterCellIndex,
1563  "[calibration] Returns the center cell index (within its 5x5 (i=5), 7x7 (i=7), 9x9 (i=9) or 11x11 (i=11) neighbours neighbours)");
1564  REGISTER_VARIABLE("eclcaldigitMCEnergy(i, j)", getMCEnergy,
1565  "[calibration] Returns the true deposited energy of all particles of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours (1-based)");
1566  REGISTER_VARIABLE("clusterNHitsThreshold(i)", getClusterNHitsThreshold,
1567  "[calibration] Returns sum of crystal weights sum(w_i) with w_i<=1 associated to this cluster above threshold (in GeV)");
1568 
1569  VARIABLE_GROUP("ECL Calibration (based on extrapolated tracks) (cDST)");
1570  REGISTER_VARIABLE("eclcaldigitExtEnergy(i, j)", getExtECLCalDigitEnergy,
1571  "[calibration] Returns the energy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1572  REGISTER_VARIABLE("eclcaldigitExtTime(i, j)", getExtECLCalDigitTime,
1573  "[calibration] Returns the time of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1574  REGISTER_VARIABLE("eclcaldigitExtTimeResolution(i, j)", getExtECLCalDigitTimeResolution,
1575  "[calibration] Returns the time resolution (dt99) of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1576  REGISTER_VARIABLE("eclcaldigitExtTwoComponentTotalEnergy(i, j)", getExtECLCalDigitTwoComponentTotalEnergy,
1577  "[calibration] Returns the TwoComponentTotalEnergy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1578  REGISTER_VARIABLE("eclcaldigitExtTwoComponentHadronEnergy(i, j)", getExtECLCalDigitTwoComponentHadronEnergy,
1579  "[calibration] Returns the TwoComponentHadronEnergy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1580  REGISTER_VARIABLE("eclcaldigitExtTwoComponentChi2(i, j)", getExtECLCalDigitTwoComponentChi2,
1581  "[calibration] Returns the TwoComponentchi2 of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1582  REGISTER_VARIABLE("eclcaldigitExtPhi(i, j)", getExtPhi,
1583  "[calibration] Returns phi of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an extrapolated track");
1584  REGISTER_VARIABLE("eclcaldigitExtTheta(i, j)", getExtTheta,
1585  "[calibration] Returns theta of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an extrapolated track");
1586  REGISTER_VARIABLE("eclcaldigitExtPhiId(i, j)", getExtPhiId,
1587  "[calibration] Returns the phi Id of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11)) neighbours for an extrapolated track");
1588  REGISTER_VARIABLE("eclcaldigitExtThetaId(i, j)", getExtThetaId,
1589  "[calibration] Returns the theta Id of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an extrapolated track");
1590  REGISTER_VARIABLE("eclcaldigitExtCellId", getExtCellId, "[calibration] Returns the extrapolated cell id");
1591  REGISTER_VARIABLE("eclcaldigitExtCellThetaId", getExtCellThetaId, "[calibration] Returns the ext cell theta id");
1592  REGISTER_VARIABLE("eclcaldigitExtCellPhiId", getExtCellPhiId, "[calibration] Returns the ext cell phi id");
1593  REGISTER_VARIABLE("eclcaldigitExtCellCrystalTheta", getExtCellCrystalTheta, "[calibration] Returns the ext cell crystal theta");
1594  REGISTER_VARIABLE("eclcaldigitExtCellCrystalPhi", getExtCellCrystalPhi, "[calibration] Returns the ext cell crystal phi");
1595  REGISTER_VARIABLE("eclcaldigitExtCenterCellIndex(i)", getExtCenterCellIndex,
1596  "[calibration] Returns the center cell index (within its 5x5 (i=5), 7x7 (i=7), 9x9 (i=9) or 11x11 (i=11) neighbours) for an ext track");
1597 
1598  REGISTER_VARIABLE("eclcaldigitExtFrontPositionPhiOffset", getExtFrontPositionPhiOffset,
1599  "[calibration] Returns the difference in the azimuthal angle (in radians)"
1600  "between the position where the track hit the front face of the ECL and the"
1601  "center of the struck crystal projected onto the front surface.");
1602  REGISTER_VARIABLE("eclcaldigitExtFrontPositionThetaOffset", getExtFrontPositionThetaOffset,
1603  "[calibration] Returns the difference in the polar angle (in radians)"
1604  "between the position where the track hit the front face of the ECL and the"
1605  "center of the struck crystal projected onto the front surface.");
1606  REGISTER_VARIABLE("eclcaldigitExtFrontPositionPhiPointing", getExtFrontPositionPhiPointing,
1607  "[calibration] Returns the difference in the azimuthal angle (in radians)"
1608  "between the momentum direction when the track hit the front face of the ECL and the"
1609  "orientation of the struck crystal.");
1610  REGISTER_VARIABLE("eclcaldigitExtFrontPositionThetaPointing", getExtFrontPositionThetaPointing,
1611  "[calibration] Returns the difference in the polar angle (in radians)"
1612  "between the momentum direction when the track hit the front face of the ECL and the"
1613  "orientation of the struck crystal.");
1614 
1615  REGISTER_VARIABLE("eclcaldigitExtTwoComponentFitType(i, j)", getExtECLCalDigitTwoComponentFitType,
1616  "[calibration] Returns the TwoComponentFitType of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1617  REGISTER_VARIABLE("eclcaldigitExtTwoComponentDiodeEnergy(i, j)", getExtECLCalDigitTwoComponentDiodeEnergy,
1618  "[calibration] Returns the TwoComponentDiodeEnergy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1619  REGISTER_VARIABLE("eclcaldigitExtTwoComponentChi2Saved_PhotonHadron(i, j)", getExtECLCalDigitTwoComponentChi2Saved_PhotonHadron,
1620  "[calibration] Returns the TwoComponentChi2Saved_PhotonHadron of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1621  REGISTER_VARIABLE("eclcaldigitExtTwoComponentChi2Saved_PileUpPhoton(i, j)", getExtECLCalDigitTwoComponentChi2Saved_PileUpPhoton,
1622  "[calibration] Returns the TwoComponentChi2Saved_PileUpPhoton of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1623  REGISTER_VARIABLE("eclcaldigitExtTwoComponentChi2Saved_PhotonDiode(i, j)", getExtECLCalDigitTwoComponentChi2Saved_PhotonDiode,
1624  "[calibration] Returns the TwoComponentChi2Saved_PhotonDiode of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours for an ext track");
1625 
1626  REGISTER_VARIABLE("eclcaldigitTwoComponentFitType(i, j)", getECLCalDigitTwoComponentFitType,
1627  "[calibration] Returns the TwoComponentFitType of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1628  REGISTER_VARIABLE("eclcaldigitTwoComponentDiodeEnergy(i, j)", getECLCalDigitTwoComponentDiodeEnergy,
1629  "[calibration] Returns the TwoComponentDiodeEnergy of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1630  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2Saved_PhotonHadron(i, j)", getECLCalDigitTwoComponentChi2Saved_PhotonHadron,
1631  "[calibration] Returns the TwoComponentChi2Saved_PhotonHadron of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1632  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2Saved_PileUpPhoton(i, j)", getECLCalDigitTwoComponentChi2Saved_PileUpPhoton,
1633  "[calibration] Returns the TwoComponentChi2Saved_PileUpPhoton of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1634  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2Saved_PhotonDiode(i, j)", getECLCalDigitTwoComponentChi2Saved_PhotonDiode,
1635  "[calibration] Returns the TwoComponentChi2Saved_PhotonDiode of the i-th caldigit for 5x5 (j=5), 7x7 (j=7), 9x9 (j=9) or 11x11 (j=11) neighbours");
1636 
1637  REGISTER_VARIABLE("eclcaldigitEnergyByEnergyRank(i)", getECLCalDigitEnergyByEnergyRank,
1638  "[calibration/eclChargedPIDExpert] Returns the caldigit energy of the i-th highest energy caldigit in the cluster (i>=0)");
1639 
1640  REGISTER_VARIABLE("eclcaldigitTimeByEnergyRank(i)", getECLCalDigitTimeByEnergyRank,
1641  "[calibration] Returns the caldigit time of the i-th highest energy caldigit in the cluster (i>=0)");
1642 
1643  REGISTER_VARIABLE("eclcaldigitTwoComponentFitTypeByEnergyRank(i[, b])", getTwoComponentFitTypeByEnergyRank,
1644  "[calibration/eclChargedPIDExpert] Returns the offline fit type of the i-th highest energy caldigit in the cluster (i>=0). If b is set to 1.0 only caldigits considered to have good quality PSD fits return PSD information.");
1645 
1646  REGISTER_VARIABLE("eclcaldigitMCEnergyByEnergyRank(i)", getECLCalDigitMCEnergyByEnergyRank,
1647  "[calibration] Returns the caldigit MC Energy of the i-th highest energy caldigit in the cluster (i>=0)");
1648 
1649  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2ByEnergyRank(i[, b])", getTwoComponentChi2ByEnergyRank,
1650  "[calibration/eclChargedPIDExpert] Returns the two component chi2 of the i-th highest energy caldigit in the cluster (i>=0). If b is set to 1.0 only caldigits considered to have good quality PSD fits return PSD information.");
1651 
1652  REGISTER_VARIABLE("eclcaldigitTwoComponentEnergyByEnergyRank(i)", getTwoComponentTotalEnergyByEnergyRank,
1653  "[calibration] Returns the two component total energy of the i-th highest energy caldigit in the cluster (i>=0)");
1654 
1655  REGISTER_VARIABLE("eclcaldigitTwoComponentHadronEnergyByEnergyRank(i[, b])", getTwoComponentHadronEnergyByEnergyRank,
1656  "[calibration/eclChargedPIDExpert] Returns the two component fit Hadron Energy of the i-th highest energy caldigit in the cluster (i>=0). If b is set to 1.0 only caldigits considered to have good quality PSD fits return PSD information.");
1657 
1658  REGISTER_VARIABLE("eclcaldigitTwoComponentDiodeEnergyByEnergyRank(i[, b])", getTwoComponentDiodeEnergyByEnergyRank,
1659  "[calibration/eclChargedPIDExpert] Returns the two component fit Diode Energy of the i-th highest energy caldigit in the cluster (i>=0). If b is set to 1.0 only caldigits considered to have good quality PSD fits return PSD information.");
1660 
1661  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2SavedByEnergyRank_PhotonHadron(i)", getTwoComponentChi2SavedByEnergyRank_PhotonHadron,
1662  "[calibration] Returns the chi2 for the photo+hadron fit type of the i-th highest energy caldigit in the cluster (i>=0)");
1663 
1664  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2SavedByEnergyRank_PileUpPhoton(i)", getTwoComponentChi2SavedByEnergyRank_PileUpPhoton,
1665  "[calibration] Returns the chi2 for the photo+hadron+pile-up photon fit type of the i-th highest energy caldigit in the cluster (i>=0)");
1666 
1667  REGISTER_VARIABLE("eclcaldigitTwoComponentChi2SavedByEnergyRank_PhotonDiode(i)", getTwoComponentChi2SavedByEnergyRank_PhotonDiode,
1668  "[calibration] Returns the chi2 for the photo+diode fit type of the i-th highest energy caldigit in the cluster (i>=0)");
1669 
1670  REGISTER_VARIABLE("eclcaldigitWeightByEnergyRank(i)", getWeightByEnergyRank,
1671  "[calibration/eclChargedPIDExpert] Returns the weight of the i-th highest energy caldigit in the cluster (i>=0)");
1672 
1673  REGISTER_VARIABLE("eclcaldigitCellIdByEnergyRank(i)", getCellIdByEnergyRank,
1674  "[calibration] Returns the cell id of the i-th highest energy caldigit in the cluster (i>=0)");
1675 
1676  REGISTER_VARIABLE("totalECLCalDigitMCEnergy", getTotalECLCalDigitMCEnergy,
1677  "[calibration] Returns total deposited MC energy in all ECLCalDigits for the MC particle");
1678 
1679  REGISTER_VARIABLE("clusterTotalECLCalDigitMCEnergy", getClusterTotalECLCalDigitMCEnergy,
1680  "[calibration] Returns total MC energy deposited in ECLCalDigits in the cluster by the related MC particle");
1681 
1682  REGISTER_VARIABLE("clusterECLCalDigitMCEnergy", getClusterECLCalDigitMCEnergy,
1683  "[calibration] Returns total deposited MC energy in all ECLCalDigits for the MC particle that are used to calculate the cluster energy");
1684 
1685  REGISTER_VARIABLE("eclcaldigitPhiByEnergyRank(i)", getPhiByEnergyRank,
1686  "Returns phi of the i-th highest energy caldigit in the cluster (i>=0)");
1687 
1688  REGISTER_VARIABLE("eclcaldigitThetaByEnergyRank(i)", getThetaByEnergyRank,
1689  "Returns theta of the i-th highest energy caldigit in the cluster (i>=0)");
1690 
1691  REGISTER_VARIABLE("eclcaldigitRByEnergyRank(i)", getRByEnergyRank,
1692  "Returns R of the i-th highest energy caldigit in the cluster (i>=0)");
1693 
1694  REGISTER_VARIABLE("eclcaldigitTwoComponentHadronEnergyFractionByEnergyRank(i[, b])",
1695  getTwoComponentHadronEnergyFractionByEnergyRank,
1696  "[eclChargedPIDExpert] Returns the hadron energy fraction of the i-th highest energy caldigit in the cluster (i>=0). If b is set to 1.0 only caldigits considered to have good quality PSD fits return PSD information.");
1697  REGISTER_VARIABLE("eclcaldigitFractionOfTotalShowerEnergyByEnergyRank(i)", getECLCalDigitFractionOfShowerEnergyByEnergyRank,
1698  "[eclChargedPIDExpert] Returns the fraction of the total Shower energy in the i-th highest energy caldigit in the Shower (i>=0). Assumes a photon hypothesis.");
1699  REGISTER_VARIABLE("eclcaldigitPhiRelativeToShowerByEnergyRank(i)", getECLCalDigitPhiRelativeToShowerByEnergyRank,
1700  "[eclChargedPIDExpert] Returns phi of the vector joining the i-th highest energy caldigit in the Shower (i>=0) to the Shower center.");
1701  REGISTER_VARIABLE("eclcaldigitThetaRelativeToShowerByEnergyRank(i)", getECLCalDigitThetaRelativeToShowerByEnergyRank,
1702  "[eclChargedPIDExpert] Returns theta of the vector joining the i-th highest energy caldigit in the Shower (i>=0) to the Shower center.");
1703  REGISTER_VARIABLE("eclcaldigitCosThetaRelativeToShowerByEnergyRank(i)", getECLCalDigitCosThetaRelativeToShowerByEnergyRank,
1704  "[eclChargedPIDExpert] Returns cos(theta) of the vector joining the i-th highest energy caldigit in the Shower (i>=0) to the Shower center.");
1705  REGISTER_VARIABLE("eclcaldigitRadiusRelativeToShowerByEnergyRank(i)", getECLCalDigitRadiusRelativeToShowerByEnergyRank,
1706  "[eclChargedPIDExpert] Returns the magnitude of the vector joining the i-th highest energy caldigit in the Shower (i>=0) to the Shower center.");
1707  }
1708 // // Create an empty module which allows basf2 to easily find the library and load it from the steering file
1709 // class EnableECLCalDigitVariablesModule: public Module {}; // Register this module to create a .map lookup file.
1710 // REG_MODULE(EnableECLCalDigitVariables); /**< register the empty module */
1712 }
void SetMagThetaPhi(DataType mag, DataType theta, DataType phi)
setter with mag, theta, phi
Definition: B2Vector3.h:259
DataType Mag() const
The magnitude (rho in spherical coordinate system).
Definition: B2Vector3.h:159
The Class for ECL Geometry Parameters.
static ECLGeometryPar * Instance()
Static method to get a reference to the ECLGeometryPar instance.
B2Vector3< double > B2Vector3D
typedef for common usage with double
Definition: B2Vector3.h:516
Abstract base class for different kinds of events.