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