Belle II Software development
SkimSampleCalculator.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#include <hlt/softwaretrigger/calculations/SkimSampleCalculator.h>
10
11#include <hlt/softwaretrigger/calculations/utilities.h>
12
13#include <analysis/utility/PCmsLabTransform.h>
14#include <analysis/utility/ReferenceFrame.h>
15
16#include <analysis/ClusterUtility/ClusterUtils.h>
17#include <analysis/ContinuumSuppression/FoxWolfram.h>
18#include <analysis/dataobjects/Particle.h>
19
20#include <analysis/variables/AcceptanceVariables.h>
21#include <analysis/variables/BelleVariables.h>
22#include <analysis/variables/ECLVariables.h>
23#include <analysis/variables/FlightInfoVariables.h>
24
25#include <mdst/dataobjects/HitPatternCDC.h>
26#include <mdst/dataobjects/KLMCluster.h>
27#include <mdst/dataobjects/PIDLikelihood.h>
28#include <mdst/dataobjects/SoftwareTriggerResult.h>
29#include <mdst/dataobjects/Track.h>
30#include <mdst/dataobjects/TrackFitResult.h>
31
32#include <cdc/dataobjects/CDCDedxTrack.h>
33
34#include <TDatabasePDG.h>
35#include <Math/Vector3D.h>
36#include <Math/Vector4D.h>
37#include <Math/VectorUtil.h>
38
39#include <cmath>
40#include <numeric>
41
42
43using namespace Belle2;
44using namespace SoftwareTrigger;
45
47 m_pionParticles("pi+:skim"), m_gammaParticles("gamma:skim"), m_pionHadParticles("pi+:hadb"), m_pionTauParticles("pi+:tau"),
48 m_KsParticles("K_S0:merged"), m_LambdaParticles("Lambda0:merged"), m_DstParticles("D*+:d0pi"), m_offIpParticles("pi+:offip"),
49 m_filterL1TrgNN("software_trigger_cut&filter&L1_trigger_nn_info"),
50 m_BpParticles("B+:BtoCharmForHLT"), m_BzParticles("B0:BtoCharmForHLT")
51{
52
53}
54
56{
57 m_pionParticles.isRequired();
58 m_gammaParticles.isRequired();
59 m_pionHadParticles.isRequired();
60 m_pionTauParticles.isRequired();
61 m_KsParticles.isOptional();
62 m_LambdaParticles.isOptional();
63 m_DstParticles.isOptional();
64 m_offIpParticles.isRequired();
65 m_BpParticles.isOptional();
66 m_BzParticles.isOptional();
67
68};
69
70void SkimSampleCalculator::doCalculation(SoftwareTriggerObject& calculationResult)
71{
72 // Prefetch some later needed objects/values
73 const Particle* gammaWithMaximumRho = getElementWithMaximumRho<Particle>(m_gammaParticles);
74 const Particle* gammaWithSecondMaximumRho = getElementWithMaximumRhoBelow<Particle>(m_gammaParticles,
75 getRho(gammaWithMaximumRho));
76 const Particle* trackWithMaximumRho = getElementWithMaximumRho<Particle>(m_pionParticles);
77 const Particle* trackWithSecondMaximumRho = getElementWithMaximumRhoBelow<Particle>(m_pionParticles,
78 getRho(trackWithMaximumRho));
79
80 const double& rhoOfECLClusterWithMaximumRho = getRhoOfECLClusterWithMaximumRho(m_pionParticles, m_gammaParticles);
81 const double& rhoOfECLClusterWithSecondMaximumRho = getRhoOfECLClusterWithMaximumRhoBelow(m_pionParticles,
83 rhoOfECLClusterWithMaximumRho);
84
85 const double& rhoOfTrackWithMaximumRho = getRho(trackWithMaximumRho);
86 const double& rhoOfTrackWithSecondMaximumRho = getRho(trackWithSecondMaximumRho);
87 const double& rhoOfGammaWithMaximumRho = getRho(gammaWithMaximumRho);
88 const double& rhoOfGammaWithSecondMaximumRho = getRho(gammaWithSecondMaximumRho);
89
90 // Simple to calculate variables
91 // EC1CMSLE
92 calculationResult["EC1CMSLE"] = rhoOfECLClusterWithMaximumRho;
93
94 // EC2CMSLE
95 calculationResult["EC2CMSLE"] = rhoOfECLClusterWithSecondMaximumRho;
96
97 // EC12CMSLE
98 calculationResult["EC12CMSLE"] = rhoOfECLClusterWithMaximumRho + rhoOfECLClusterWithSecondMaximumRho;
99
100 // nTracksLE
101 calculationResult["nTracksLE"] = m_pionParticles->getListSize();
102
103 // nTracksTAU
104 calculationResult["nTracksTAU"] = m_pionTauParticles->getListSize();
105
106 // nGammasLE
107 calculationResult["nGammasLE"] = m_gammaParticles->getListSize();
108
109 // P1CMSBhabhaLE
110 calculationResult["P1CMSBhabhaLE"] = rhoOfTrackWithMaximumRho;
111
112 // P1CMSBhabhaLE/E_beam
113 calculationResult["P1OEbeamCMSBhabhaLE"] = rhoOfTrackWithMaximumRho / BeamEnergyCMS();
114
115 // P2CMSBhabhaLE
116 calculationResult["P2CMSBhabhaLE"] = rhoOfTrackWithSecondMaximumRho;
117
118 // P2CMSBhabhaLE/E_beam
119 calculationResult["P2OEbeamCMSBhabhaLE"] = rhoOfTrackWithSecondMaximumRho / BeamEnergyCMS();
120
121 // P12CMSBhabhaLE
122 calculationResult["P12CMSBhabhaLE"] = rhoOfTrackWithMaximumRho + rhoOfTrackWithSecondMaximumRho;
123
124 //G1CMSLE, the largest energy of gamma in CMS
125 calculationResult["G1CMSBhabhaLE"] = rhoOfGammaWithMaximumRho;
126 //G1OEbeamCMSLE, the largest energy of gamma in CMS over beam energy
127 calculationResult["G1OEbeamCMSBhabhaLE"] = rhoOfGammaWithMaximumRho / BeamEnergyCMS();
128
129 //G2CMSLE, the secondary largest energy of gamma in CMS
130 calculationResult["G2CMSBhabhaLE"] = rhoOfGammaWithSecondMaximumRho;
131 //G2OEbeamCMSLE, the largest energy of gamma in CMS over beam energy
132 calculationResult["G2OEbeamCMSBhabhaLE"] = rhoOfGammaWithSecondMaximumRho / BeamEnergyCMS();
133
134 //G12CMSLE, the secondary largest energy of gamma in CMS
135 calculationResult["G12CMSBhabhaLE"] = rhoOfGammaWithMaximumRho + rhoOfGammaWithSecondMaximumRho;
136 //G12CMSLE, the secondary largest energy of gamma in CMS over beam energy
137 calculationResult["G12OEbeamCMSBhabhaLE"] =
138 (rhoOfGammaWithMaximumRho + rhoOfGammaWithSecondMaximumRho) / BeamEnergyCMS();
139
140
141 // Medium hard to calculate variables
142 // ENeutralLE
143 if (gammaWithMaximumRho) {
144 calculationResult["ENeutralLE"] = getRho(gammaWithMaximumRho);
145 } else {
146 calculationResult["ENeutralLE"] = -1;
147 }
148
149 // nECLMatchTracksLE
150 const unsigned int numberOfTracksWithECLMatch = std::count_if(m_pionParticles->begin(), m_pionParticles->end(),
151 [](const Particle & particle) {
152 return particle.getECLCluster() != nullptr;
153 });
154 calculationResult["nECLMatchTracksLE"] = numberOfTracksWithECLMatch;
155
156 //nECLClustersLE
157 double neclClusters = -1.;
158 double eneclClusters = 0.;
159 StoreArray<ECLCluster> eclClusters;
160 ClusterUtils Cl;
161 double PzGamma = 0.;
162 double EsumGamma = 0.;
163 if (eclClusters.isValid()) {
164 const unsigned int numberOfECLClusters = std::count_if(eclClusters.begin(), eclClusters.end(),
165 [](const ECLCluster & eclcluster) {
166 return (eclcluster.hasHypothesis(
167 ECLCluster::EHypothesisBit::c_nPhotons)
168 and eclcluster.getEnergy(
169 ECLCluster::EHypothesisBit::c_nPhotons) > 0.1);
170 });
171 neclClusters = numberOfECLClusters;
172
173 for (int ncl = 0; ncl < eclClusters.getEntries(); ncl++) {
174 if (eclClusters[ncl]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)
175 && eclClusters[ncl]->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) > 0.1) {
176 eneclClusters += eclClusters[ncl]->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons);
177 if (!eclClusters[ncl]->getRelatedFrom<Track>()) {
178 ROOT::Math::PxPyPzEVector V4Gamma_CMS = PCmsLabTransform::labToCms(Cl.Get4MomentumFromCluster(eclClusters[ncl],
180 EsumGamma += V4Gamma_CMS.E();
181 PzGamma += V4Gamma_CMS.Pz();
182 }
183 }
184 }
185 }
186 calculationResult["nECLClustersLE"] = neclClusters;
187
188 int nb2bcc_PhiHigh = 0;
189 int nb2bcc_PhiLow = 0;
190 int nb2bcc_3D = 0;
191 ClusterUtils C;
192 for (int i = 0; i < eclClusters.getEntries() - 1; i++) {
193 if (!eclClusters[i]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
194 continue;
195 ROOT::Math::PxPyPzEVector V4g1 = C.Get4MomentumFromCluster(eclClusters[i], ECLCluster::EHypothesisBit::c_nPhotons);
196 double Eg1 = V4g1.E();
197 for (int j = i + 1; j < eclClusters.getEntries(); j++) {
198 if (!eclClusters[j]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
199 continue;
200 ROOT::Math::PxPyPzEVector V4g2 = C.Get4MomentumFromCluster(eclClusters[j], ECLCluster::EHypothesisBit::c_nPhotons);
201 double Eg2 = V4g2.E();
202 const ROOT::Math::PxPyPzEVector V4g1CMS = PCmsLabTransform::labToCms(V4g1);
203 const ROOT::Math::PxPyPzEVector V4g2CMS = PCmsLabTransform::labToCms(V4g2);
204 double Thetag1 = V4g1CMS.Theta() * TMath::RadToDeg();
205 double Thetag2 = V4g2CMS.Theta() * TMath::RadToDeg();
206 double deltphi = std::abs(ROOT::Math::VectorUtil::DeltaPhi(V4g1CMS, V4g2CMS) * TMath::RadToDeg());
207 double Tsum = Thetag1 + Thetag2;
208 if (deltphi > 170. && (Eg1 > 0.25 && Eg2 > 0.25)) nb2bcc_PhiHigh++;
209 if (deltphi > 170. && (Eg1 < 0.25 || Eg2 < 0.25)) nb2bcc_PhiLow++;
210 if (deltphi > 160. && (Tsum > 160. && Tsum < 200.)) nb2bcc_3D++;
211 }
212 }
213
214 calculationResult["nB2BCCPhiHighLE"] = nb2bcc_PhiHigh;
215 calculationResult["nB2BCCPhiLowLE"] = nb2bcc_PhiLow;
216 calculationResult["nB2BCC3DLE"] = nb2bcc_3D;
217
218
219 // AngleGTLE
220 double angleGTLE = -10.;
221 if (gammaWithMaximumRho) {
222 const ROOT::Math::XYZVector& V3g1 = gammaWithMaximumRho->getMomentum();
223 if (trackWithMaximumRho) {
224 const ROOT::Math::XYZVector& V3p1 = trackWithMaximumRho->getMomentum();
225 const double theta1 = ROOT::Math::VectorUtil::Angle(V3g1, V3p1);
226 if (angleGTLE < theta1) angleGTLE = theta1;
227 }
228 if (trackWithSecondMaximumRho) {
229 const ROOT::Math::XYZVector& V3p2 = trackWithSecondMaximumRho->getMomentum();
230 const double theta2 = ROOT::Math::VectorUtil::Angle(V3g1, V3p2);
231 if (angleGTLE < theta2) angleGTLE = theta2;
232 }
233 }
234
235 calculationResult["AngleGTLE"] = angleGTLE;
236
237 // AngleG1G2LE
238 double angleG1G2CMSLE = -10.;
239 if (gammaWithMaximumRho) {
240 const ROOT::Math::PxPyPzEVector& V4p1 = PCmsLabTransform::labToCms(gammaWithMaximumRho->get4Vector());
241 if (gammaWithSecondMaximumRho) {
242 const ROOT::Math::PxPyPzEVector& V4p2 = PCmsLabTransform::labToCms(gammaWithSecondMaximumRho->get4Vector());
243 angleG1G2CMSLE = ROOT::Math::VectorUtil::Angle(V4p1, V4p2);
244 }
245 }
246
247 calculationResult["AngleG1G2CMSLE"] = angleG1G2CMSLE;
248
249 // maxAngleTTLE
250 double maxAngleTTLE = -10.;
251 int nJpsi = 0;
252 double Jpsi = 0.;
253 const double jPsiMasswindow = 0.11;
254 if (m_pionParticles->getListSize() >= 2) {
255 for (unsigned int i = 0; i < m_pionParticles->getListSize() - 1; i++) {
256 Particle* par1 = m_pionParticles->getParticle(i);
257 for (unsigned int j = i + 1; j < m_pionParticles->getListSize(); j++) {
258 Particle* par2 = m_pionParticles->getParticle(j);
259 ROOT::Math::PxPyPzEVector V4p1 = par1->get4Vector();
260 ROOT::Math::PxPyPzEVector V4p2 = par2->get4Vector();
261 ROOT::Math::PxPyPzEVector V4pSum = V4p1 + V4p2;
262 const auto chSum = par1->getCharge() + par2->getCharge();
263 const double mSum = V4pSum.M();
264 const double JpsidM = mSum - TDatabasePDG::Instance()->GetParticle(443)->Mass();
265 if (std::abs(JpsidM) < jPsiMasswindow && chSum == 0) nJpsi++;
266 const ROOT::Math::PxPyPzEVector V4p1CMS = PCmsLabTransform::labToCms(V4p1);
267 const ROOT::Math::PxPyPzEVector V4p2CMS = PCmsLabTransform::labToCms(V4p2);
268 const double temp = ROOT::Math::VectorUtil::Angle(V4p1CMS, V4p2CMS);
269 if (maxAngleTTLE < temp) maxAngleTTLE = temp;
270 }
271 }
272 }
273
274 if (nJpsi != 0) Jpsi = 1;
275
276 calculationResult["maxAngleTTLE"] = maxAngleTTLE;
277 calculationResult["Jpsi"] = Jpsi;
278
279 //maxAngleGGLE
280 double maxAngleGGLE = -10.;
281 if (m_gammaParticles->getListSize() >= 2) {
282 for (unsigned int i = 0; i < m_gammaParticles->getListSize() - 1; i++) {
283 Particle* par1 = m_gammaParticles->getParticle(i);
284 for (unsigned int j = i + 1; j < m_gammaParticles->getListSize(); j++) {
285 Particle* par2 = m_gammaParticles->getParticle(j);
286 ROOT::Math::PxPyPzEVector V4p1 = PCmsLabTransform::labToCms(par1->get4Vector());
287 ROOT::Math::PxPyPzEVector V4p2 = PCmsLabTransform::labToCms(par2->get4Vector());
288 const double temp = ROOT::Math::VectorUtil::Angle(V4p1, V4p2);
289 if (maxAngleGGLE < temp) maxAngleGGLE = temp;
290 }
291 }
292 }
293
294 calculationResult["maxAngleGGLE"] = maxAngleGGLE;
295
296 // nEidLE
297 const unsigned int nEidLE = std::count_if(m_pionParticles->begin(), m_pionParticles->end(),
298 [](const Particle & p) {
299 const double& momentum = p.getMomentumMagnitude();
300 const double& r_rho = getRho(&p);
301 const ECLCluster* eclTrack = p.getECLCluster();
302 if (eclTrack) {
303 const double& energyOverMomentum = eclTrack->getEnergy(
304 ECLCluster::EHypothesisBit::c_nPhotons) / momentum;
305 double r_rhotoebeam = r_rho / BeamEnergyCMS();
306 return (r_rhotoebeam) > 0.35 && energyOverMomentum > 0.8;
307 }
308 return false;
309 });
310
311 calculationResult["nEidLE"] = nEidLE;
312
313
314 // VisibleEnergyLE
315 const double visibleEnergyTracks = std::accumulate(m_pionParticles->begin(), m_pionParticles->end(), 0.0,
316 [](const double & visibleEnergy, const Particle & p) {
317 return visibleEnergy + p.getMomentumMagnitude();
318 });
319
320 const double visibleEnergyGammas = std::accumulate(m_gammaParticles->begin(), m_gammaParticles->end(), 0.0,
321 [](const double & visibleEnergy, const Particle & p) {
322 return visibleEnergy + p.getMomentumMagnitude();
323 });
324
325 calculationResult["VisibleEnergyLE"] = visibleEnergyTracks + visibleEnergyGammas;
326
327 // EtotLE
328 const double eTotTracks = std::accumulate(m_pionParticles->begin(), m_pionParticles->end(), 0.0,
329 [](const double & eTot, const Particle & p) {
330 const ECLCluster* eclCluster = p.getECLCluster();
331 if (eclCluster) {
332 const double eclEnergy = eclCluster->getEnergy(
333 ECLCluster::EHypothesisBit::c_nPhotons);
334 if (eclEnergy > 0.1) {
335 return eTot + eclCluster->getEnergy(
336 ECLCluster::EHypothesisBit::c_nPhotons);
337 }
338 }
339 return eTot;
340 });
341
342 const double eTotGammas = std::accumulate(m_gammaParticles->begin(), m_gammaParticles->end(), 0.0,
343 [](const double & eTot, const Particle & p) {
344 return eTot + p.getEnergy();
345 });
346 double Etot = eTotTracks + eTotGammas;
347 calculationResult["EtotLE"] = Etot;
348
349 //KLM information
350 // The clusters with the largest pentrate layers in KLM.
351 double numMaxLayerKLM = -1.;
352 double numSecMaxLayerKLM = -1.;
353 StoreArray<KLMCluster> klmClusters;
354 if (klmClusters.isValid()) {
355 for (const auto& klmCluster : klmClusters) {
356 double klmClusterLayer = klmCluster.getLayers();
357 if (numMaxLayerKLM < klmClusterLayer) {
358 numSecMaxLayerKLM = numMaxLayerKLM;
359 numMaxLayerKLM = klmClusterLayer;
360 } else if (numSecMaxLayerKLM < klmClusterLayer)
361 numSecMaxLayerKLM = klmClusterLayer;
362 }
363 }
364 calculationResult["N1KLMLayer"] = numMaxLayerKLM;
365 calculationResult["N2KLMLayer"] = numSecMaxLayerKLM;
366
367 //define bhabha_2trk, bhabha_1trk, eclbhabha
368 int charget1 = -10;
369 if (trackWithMaximumRho) charget1 = trackWithMaximumRho->getCharge();
370 int charget2 = -10;
371 if (trackWithSecondMaximumRho) charget2 = trackWithSecondMaximumRho->getCharge();
372
373 double Bhabha2Trk = 0.;
374 int ntrk_bha = m_pionParticles->getListSize();
375 double rp1ob = rhoOfTrackWithMaximumRho / BeamEnergyCMS();
376 double rp2ob = rhoOfTrackWithSecondMaximumRho / BeamEnergyCMS();
377 bool bhabha2trk_tag =
378 ntrk_bha >= 2 && maxAngleTTLE > 2.88 && nEidLE >= 1 && rp1ob > 0.35 && rp2ob > 0.35 && (Etot) > 4.0
379 && (std::abs(charget1) == 1 && std::abs(charget2) == 1 && (charget1 + charget2) == 0);
380 if (bhabha2trk_tag) Bhabha2Trk = 1;
381 calculationResult["Bhabha2Trk"] = Bhabha2Trk;
382
383 double Bhabha1Trk = 0.;
384 double rc1ob = rhoOfGammaWithMaximumRho / BeamEnergyCMS();
385 double rc2ob = rhoOfGammaWithSecondMaximumRho / BeamEnergyCMS();
386 bool bhabha1trk_tag = ntrk_bha == 1 && rp1ob > 0.35 && rc1ob > 0.35 && angleGTLE > 2.618;
387 if (bhabha1trk_tag) Bhabha1Trk = 1;
388 calculationResult["Bhabha1Trk"] = Bhabha1Trk;
389
390 double ggSel = 0.;
391 bool gg_tag = ntrk_bha <= 1 && nEidLE == 0 && rc1ob > 0.35 && rc2ob > 0.2 && Etot > 4.0 && maxAngleGGLE > 2.618;
392 if (gg_tag) ggSel = 1;
393 calculationResult["GG"] = ggSel;
394
395 // Bhabha skim with ECL information only (bhabhaecl)
396 double BhabhaECL = 0.;
397 ClusterUtils Cls;
398 for (int i = 0; i < eclClusters.getEntries() - 1; i++) {
399 if (!eclClusters[i]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
400 continue;
401
402 ROOT::Math::PxPyPzEVector V4g1 = PCmsLabTransform::labToCms(Cls.Get4MomentumFromCluster(eclClusters[i],
404 double Eg1ob = V4g1.E() / (2 * BeamEnergyCMS());
405 for (int j = i + 1; j < eclClusters.getEntries(); j++) {
406 if (!eclClusters[j]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
407 continue;
408 ROOT::Math::PxPyPzEVector V4g2 = PCmsLabTransform::labToCms(Cls.Get4MomentumFromCluster(eclClusters[j],
410 double Eg2ob = V4g2.E() / (2 * BeamEnergyCMS());
411 double Thetag1 = V4g1.Theta() * TMath::RadToDeg();
412 double Thetag2 = V4g2.Theta() * TMath::RadToDeg();
413 double deltphi = std::abs(ROOT::Math::VectorUtil::DeltaPhi(V4g1, V4g2) * TMath::RadToDeg());
414 double Tsum = Thetag1 + Thetag2;
415 if ((deltphi > 165. && deltphi < 178.5) && (Eg1ob > 0.4 && Eg2ob > 0.4 && (Eg1ob > 0.45 || Eg2ob > 0.45)) && (Tsum > 178.
416 && Tsum < 182.)) BhabhaECL = 1;
417 }
418 }
419 calculationResult["BhabhaECL"] = BhabhaECL;
420
421 // Bhabha skim (BhabhaCDC) for CDC dE/dx calib studies
422 double BhabhaCDC = 0.;
423 // Radiative Bhabha skim (radee)
424 double radee = 0.;
425 const double lowdEdxEdge = 0.70, highdEdxEdge = 1.30;
426 const double lowEoPEdge = 0.70, highEoPEdge = 1.30;
427
428 if (m_pionParticles->getListSize() == 2) {
429
430 //------------First track variables----------------
431 for (unsigned int i = 0; i < m_pionParticles->getListSize() - 1; i++) {
432
433 Particle* part1 = m_pionParticles->getParticle(i);
434 if (!part1) continue;
435
436 const auto chargep1 = part1->getCharge();
437 if (std::abs(chargep1) != 1) continue;
438
439 const ECLCluster* eclTrack1 = part1->getECLCluster();
440 if (!eclTrack1) continue;
441 if (!eclTrack1->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)) continue;
442
443 const double& momentum1 = part1->getMomentumMagnitude();
444 const double& energyOverMomentum1 = eclTrack1->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) / momentum1;
445 if (energyOverMomentum1 <= lowEoPEdge || energyOverMomentum1 >= highEoPEdge) continue;
446
447 const Track* track1 = part1->getTrack();
448 if (!track1) continue;
449
450 const TrackFitResult* trackFit1 = track1->getTrackFitResultWithClosestMass(Const::pion);
451 if (!trackFit1) continue;
452 if (trackFit1->getHitPatternCDC().getNHits() <= 0) continue;
453
454 //------------Second track variables----------------
455 for (unsigned int j = i + 1; j < m_pionParticles->getListSize(); j++) {
456
457 Particle* part2 = m_pionParticles->getParticle(j);
458 if (!part2) continue;
459
460 const auto chargep2 = part2->getCharge();
461 if (std::abs(chargep2) != 1 || (chargep1 + chargep2 != 0)) continue;
462
463 const ECLCluster* eclTrack2 = part2->getECLCluster();
464 if (!eclTrack2) continue;
465 if (!eclTrack2->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)) continue;
466
467 const double& momentum2 = part2->getMomentumMagnitude();
468 const double& energyOverMomentum2 = eclTrack2->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) / momentum2;
469 if (energyOverMomentum2 <= lowEoPEdge || energyOverMomentum2 >= highEoPEdge) continue;
470
471 const Track* track2 = part2->getTrack();
472 if (!track2) continue;
473
474 const TrackFitResult* trackFit2 = track2->getTrackFitResultWithClosestMass(Const::pion);
475 if (!trackFit2) continue;
476 if (trackFit2->getHitPatternCDC().getNHits() <= 0) continue;
477
478 BhabhaCDC = 1;
479
480 const CDCDedxTrack* dedxTrack1 = track1->getRelatedTo<CDCDedxTrack>();
481 if (!dedxTrack1) continue;
482
483 const CDCDedxTrack* dedxTrack2 = track2->getRelatedTo<CDCDedxTrack>();
484 if (!dedxTrack2) continue;
485
486 double p1_dedxnosat = dedxTrack1->getDedxNoSat();
487 double p2_dedxnosat = dedxTrack2->getDedxNoSat();
488
489 if ((p1_dedxnosat > lowdEdxEdge && p1_dedxnosat < highdEdxEdge) || (p2_dedxnosat > lowdEdxEdge
490 && p2_dedxnosat < highdEdxEdge)) radee = 1;
491
492 }
493 }
494 }
495
496 calculationResult["BhabhaCDC"] = BhabhaCDC;
497 calculationResult["Radee"] = radee;
498
499 // Dimuon skim (mumutight) taken from the offline skim + Radiative dimuon (radmumu)
500 double mumutight = 0.;
501 double eMumuTotGammas = 0.;
502 int nTracks = 0;
503 double radmumu = 0.;
504 int nGammas = m_gammaParticles->getListSize();
505
506 for (int t = 0; t < nGammas; t++) {
507 const Particle* part = m_gammaParticles->getParticle(t);
508 const auto& frame = ReferenceFrame::GetCurrent();
509 eMumuTotGammas += frame.getMomentum(part).E();
510 }
511
512 StoreArray<Track> tracks;
513 nTracks = tracks.getEntries();
514 PCmsLabTransform T;
515 const ROOT::Math::PxPyPzEVector pIN = T.getBeamFourMomentum();
516 const auto& fr = ReferenceFrame::GetCurrent();
517
518 if (m_pionParticles->getListSize() == 2) {
519
520 //------------First track variables----------------
521 for (unsigned int k = 0; k < m_pionParticles->getListSize() - 1; k++) {
522
523 Particle* part1 = m_pionParticles->getParticle(k);
524 if (!part1) continue;
525
526 const auto chargep1 = part1->getCharge();
527 if (std::abs(chargep1) != 1) continue;
528
529 const ECLCluster* eclTrack1 = part1->getECLCluster();
530 if (!eclTrack1) continue;
531 if (!eclTrack1->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)) continue;
532
533 const Track* track1 = part1->getTrack();
534 if (!track1) continue;
535
536 const TrackFitResult* trackFit1 = track1->getTrackFitResultWithClosestMass(Const::pion);
537 if (!trackFit1) continue;
538
539 const ROOT::Math::PxPyPzEVector V4p1 = trackFit1->get4Momentum();
540 const ROOT::Math::PxPyPzEVector V4p1CMS = PCmsLabTransform::labToCms(V4p1);
541
542 const double p1MomLab = V4p1.P();
543 double highestP = p1MomLab;
544 const double p1CDChits = trackFit1->getHitPatternCDC().getNHits();
545 const PIDLikelihood* p1Pid = part1->getPIDLikelihood();
546 bool p1hasKLMid = 0;
547 if (p1Pid) p1hasKLMid = p1Pid->isAvailable(Const::KLM);
548 const double p1isInCDC = Variable::inCDCAcceptance(part1);
549 const double p1clusPhi = Variable::eclClusterPhi(part1);
550
551 const double Pp1 = V4p1CMS.R();
552 const double Thetap1 = V4p1CMS.Theta() * TMath::RadToDeg();
553 const double Phip1 = V4p1CMS.Phi() * TMath::RadToDeg();
554
555 const double enECLTrack1 = eclTrack1->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons);
556
557 const bool goodTrk1 = enECLTrack1 > 0 && enECLTrack1 < 0.5 && p1CDChits > 0
558 && ((p1hasKLMid == 0 && enECLTrack1 < 0.25 && p1MomLab < 2.0) || p1hasKLMid == 1) && p1isInCDC == 1;
559
560 //------------Second track variables----------------
561 for (unsigned int l = k + 1; l < m_pionParticles->getListSize(); l++) {
562
563 Particle* part2 = m_pionParticles->getParticle(l);
564 if (!part2) continue;
565
566 const auto chargep2 = part2->getCharge();
567 if (std::abs(chargep2) != 1 || (chargep1 + chargep2 != 0)) continue;
568
569 const ECLCluster* eclTrack2 = part2->getECLCluster();
570 if (!eclTrack2) continue;
571 if (!eclTrack2->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)) continue;
572
573 const Track* track2 = part2->getTrack();
574 if (!track2) continue;
575
576 const TrackFitResult* trackFit2 = track2->getTrackFitResultWithClosestMass(Const::pion);
577 if (!trackFit2) continue;
578
579 const ROOT::Math::PxPyPzEVector V4p2 = trackFit2->get4Momentum();
580 const ROOT::Math::PxPyPzEVector V4p2CMS = PCmsLabTransform::labToCms(V4p2);
581
582 const double p2MomLab = V4p2.P();
583 double lowestP = p2MomLab;
584 const double p2CDChits = trackFit2->getHitPatternCDC().getNHits();
585 const PIDLikelihood* p2Pid = part2->getPIDLikelihood();
586 bool p2hasKLMid = 0;
587 if (p2Pid) p2hasKLMid = p2Pid->isAvailable(Const::KLM);
588 const double p2isInCDC = Variable::inCDCAcceptance(part2);
589 const double p2clusPhi = Variable::eclClusterPhi(part2);
590
591 const double Pp2 = V4p2CMS.R();
592 const double Thetap2 = V4p2CMS.Theta() * TMath::RadToDeg();
593 const double Phip2 = V4p2CMS.Phi() * TMath::RadToDeg();
594
595 const double acopPhi = std::abs(180 - std::abs(Phip1 - Phip2));
596 const double acopTheta = std::abs(std::abs(Thetap1 + Thetap2) - 180);
597
598 const double enECLTrack2 = eclTrack2->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons);
599
600 const bool goodTrk2 = enECLTrack2 > 0 && enECLTrack2 < 0.5 && p2CDChits > 0
601 && ((p2hasKLMid == 0 && enECLTrack2 < 0.25 && p2MomLab < 2.0) || p2hasKLMid == 1) && p2isInCDC == 1;
602
603 double eTotMumuTracks = enECLTrack1 + enECLTrack2;
604 double EMumutot = eTotMumuTracks + eMumuTotGammas;
605
606 bool mumutight_tag = enECLTrack1 < 0.5 && enECLTrack2 < 0.5 && EMumutot < 2 && acopPhi < 10 && acopTheta < 10 && nTracks == 2
607 && Pp1 > 0.5 && Pp2 > 0.5;
608
609 if (mumutight_tag) mumutight = 1;
610
611 if (p1MomLab < p2MomLab) {
612 lowestP = highestP;
613 highestP = p2MomLab;
614 }
615
616 double diffPhi = p1clusPhi - p2clusPhi;
617 if (std::abs(diffPhi) > M_PI) {
618 if (diffPhi > M_PI) {
619 diffPhi = diffPhi - 2 * M_PI;
620 } else {
621 diffPhi = 2 * M_PI + diffPhi;
622 }
623 }
624
625 const double recoilP = fr.getMomentum(pIN - V4p1 - V4p2).P();
626
627 const bool radmumu_tag = nTracks < 4 && goodTrk1 == 1 && goodTrk2 == 1 && highestP > 1 && lowestP < 3
628 && (p1hasKLMid == 1 || p2hasKLMid == 1) && std::abs(diffPhi) >= 0.5 * M_PI && recoilP > 0.1
629 && (enECLTrack1 <= 0.25 || enECLTrack2 <= 0.25);
630
631 if (radmumu_tag) radmumu = 1;
632
633 }
634 }
635 }
636
637 calculationResult["MumuTight"] = mumutight;
638 calculationResult["Radmumu"] = radmumu;
639
640 //Retrieve variables for HadronB skims
641 double EsumPiHad = 0;
642 double PzPiHad = 0;
643 int nHadTracks = m_pionHadParticles->getListSize();
644 double hadronb = 0;
645 double hadronb1 = 0;
646 double hadronb2 = 0;
647 std::vector<ROOT::Math::PxPyPzEVector> m_pionHad;
648
649 for (int nPiHad = 0; nPiHad < nHadTracks; nPiHad++) {
650 Particle* parPiHad = m_pionHadParticles->getParticle(nPiHad);
651 ROOT::Math::PxPyPzEVector V4PiHad = PCmsLabTransform::labToCms(parPiHad->get4Vector());
652 m_pionHad.push_back(V4PiHad);
653 EsumPiHad += V4PiHad.E();
654 PzPiHad += V4PiHad.Pz();
655 }
656
657 double visibleEnergyCMSnorm = (EsumPiHad + EsumGamma) / (BeamEnergyCMS() * 2.0);
658 double EsumCMSnorm = eneclClusters / (BeamEnergyCMS() * 2.0);
659 double PzTotCMSnorm = (PzPiHad + PzGamma) / (BeamEnergyCMS() * 2.0);
660
661 bool hadronb_tag = nHadTracks >= 3 && visibleEnergyCMSnorm > 0.2 && std::abs(PzTotCMSnorm) < 0.5 && neclClusters > 1
662 && EsumCMSnorm > 0.1 && EsumCMSnorm < 0.8;
663
664 if (hadronb_tag) {
665 hadronb = 1;
666 FoxWolfram fw(m_pionHad);
667 fw.calculateBasicMoments();
668 double R2 = fw.getR(2);
669 if (R2 < 0.4) hadronb1 = 1;
670 if (hadronb1 && nHadTracks >= 5) hadronb2 = 1;
671 }
672
673 calculationResult["HadronB"] = hadronb;
674 calculationResult["HadronB1"] = hadronb1;
675 calculationResult["HadronB2"] = hadronb2;
676
677 // nKshort
678 int nKshort = 0;
679 double Kshort = 0.;
680 const double KsMassLow = 0.468;
681 const double KsMassHigh = 0.528;
682
683 if (m_KsParticles.isValid()) {
684 for (unsigned int i = 0; i < m_KsParticles->getListSize(); i++) {
685 const Particle* mergeKsCand = m_KsParticles->getParticle(i);
686 const double isKsCandGood = Variable::goodBelleKshort(mergeKsCand);
687 const double KsCandMass = mergeKsCand->getMass();
688 if (KsCandMass > KsMassLow && KsCandMass < KsMassHigh && isKsCandGood == 1.) nKshort++;
689 }
690 }
691
692 if (nKshort != 0) Kshort = 1;
693
694 calculationResult["Kshort"] = Kshort;
695
696 // 4 leptons skim
697 int nFourLep = 0;
698 double fourLep = 0.;
699
700 const double visibleEnergyCMS = visibleEnergyCMSnorm * BeamEnergyCMS() * 2.0;
701 const unsigned int n_particles = m_pionHadParticles->getListSize();
702
703 if (n_particles >= 2) {
704 for (unsigned int i = 0; i < n_particles - 1; i++) {
705 Particle* par1 = m_pionHadParticles->getParticle(i);
706 for (unsigned int j = i + 1; j < n_particles; j++) {
707 Particle* par2 = m_pionHadParticles->getParticle(j);
708 const auto chSum = par1->getCharge() + par2->getCharge();
709 const ROOT::Math::PxPyPzEVector V4p1 = par1->get4Vector();
710 const ROOT::Math::PxPyPzEVector V4p2 = par2->get4Vector();
711 const double opAng = V4p1.Theta() + V4p2.Theta();
712 const ROOT::Math::PxPyPzEVector V4pSum = V4p1 + V4p2;
713 const ROOT::Math::PxPyPzEVector V4pSumCMS = PCmsLabTransform::labToCms(V4pSum);
714 const double ptCMS = V4pSumCMS.Pt();
715 const double pzCMS = V4pSumCMS.Pz();
716 const double mSum = V4pSum.M();
717
718 const bool fourLepCand = chSum == 0 && (V4p1.P() > 0.4 && V4p2.P() > 0.4) && cos(opAng) > -0.997 && ptCMS < 0.15
719 && std::abs(pzCMS) < 2.5
720 && mSum < 6;
721
722 if (fourLepCand) nFourLep++;
723 }
724 }
725 }
726
727 if (nFourLep != 0 && visibleEnergyCMS < 6) fourLep = 1;
728
729 calculationResult["FourLep"] = fourLep;
730
731 // nLambda
732 unsigned int nLambda = 0;
733
734 if (m_LambdaParticles.isValid()) {
735 for (unsigned int i = 0; i < m_LambdaParticles->getListSize(); i++) {
736 const Particle* mergeLambdaCand = m_LambdaParticles->getParticle(i);
737 const double flightDist = Variable::flightDistance(mergeLambdaCand);
738 const double flightDistErr = Variable::flightDistanceErr(mergeLambdaCand);
739 const double flightSign = flightDist / flightDistErr;
740 const Particle* protCand = mergeLambdaCand->getDaughter(0);
741 const Particle* pionCand = mergeLambdaCand->getDaughter(1);
742 const double protMom = protCand->getP();
743 const double pionMom = pionCand->getP();
744 const double asymPDaughters = (protMom - pionMom) / (protMom + pionMom);
745 if (flightSign > 10 && asymPDaughters > 0.41) nLambda++;
746 }
747 }
748
749 if (nLambda > 0) {
750 calculationResult["Lambda"] = 1;
751 } else {
752 calculationResult["Lambda"] = 0;
753 }
754
755 // nDstp
756 unsigned int nDstp1 = 0;
757 unsigned int nDstp2 = 0;
758 unsigned int nDstp3 = 0;
759 unsigned int nDstp4 = 0;
760
761 if (m_DstParticles.isValid() && (ntrk_bha >= 3 && Bhabha2Trk == 0)) {
762 for (unsigned int i = 0; i < m_DstParticles->getListSize(); i++) {
763 const Particle* allDstCand = m_DstParticles->getParticle(i);
764 const double dstDecID = allDstCand->getExtraInfo("decayModeID");
765 if (dstDecID == 1.) nDstp1++;
766 if (dstDecID == 2.) nDstp2++;
767 if (dstDecID == 3.) nDstp3++;
768 if (dstDecID == 4.) nDstp4++;
769 }
770 }
771
772
773 if (nDstp1 > 0) {
774 calculationResult["Dstp1"] = 1;
775 } else {
776 calculationResult["Dstp1"] = 0;
777 }
778
779 if (nDstp2 > 0) {
780 calculationResult["Dstp2"] = 1;
781 } else {
782 calculationResult["Dstp2"] = 0;
783 }
784
785 if (nDstp3 > 0) {
786 calculationResult["Dstp3"] = 1;
787 } else {
788 calculationResult["Dstp3"] = 0;
789 }
790
791 if (nDstp4 > 0) {
792 calculationResult["Dstp4"] = 1;
793 } else {
794 calculationResult["Dstp4"] = 0;
795 }
796
797 // nTracksOffIP
798 calculationResult["nTracksOffIP"] = m_offIpParticles->getListSize();
799
800 // Flag for events with Trigger B2Link information
801 calculationResult["NeuroTRG"] = 0;
802 calculationResult["GammaGammaFilter"] = 0;
803
804 StoreObjPtr<SoftwareTriggerResult> filter_result;
805 if (filter_result.isValid()) {
806 const std::map<std::string, int>& nonPrescaledResults = filter_result->getNonPrescaledResults();
807 if (nonPrescaledResults.find(m_filterL1TrgNN) != nonPrescaledResults.end()) {
808 const bool hasNN = (filter_result->getNonPrescaledResult(m_filterL1TrgNN) == SoftwareTriggerCutResult::c_accept);
809 if (hasNN) calculationResult["NeuroTRG"] = 1;
810 }
811 const bool ggEndcap = (filter_result->getNonPrescaledResult("software_trigger_cut&filter&ggEndcapLoose") ==
813 const bool ggBarrel = (filter_result->getNonPrescaledResult("software_trigger_cut&filter&ggBarrelLoose") ==
815 if (ggEndcap || ggBarrel) calculationResult["GammaGammaFilter"] = 1;
816 }
817
818 //Dimuon skim with invariant mass cut allowing at most one track not to be associated with ECL clusters
819
820 double mumuHighMass = 0.;
821
822 if (trackWithMaximumRho && trackWithSecondMaximumRho) {
823 int hasClus = 0;
824 double eclE1 = 0., eclE2 = 0.;
825
826 const auto charge1 = trackWithMaximumRho->getCharge();
827 const auto charge2 = trackWithSecondMaximumRho->getCharge();
828 const auto chSum = charge1 + charge2;
829
830 const ECLCluster* eclTrack1 = trackWithMaximumRho->getECLCluster();
831 if (eclTrack1) {
832 hasClus++;
834 }
835
836 const ECLCluster* eclTrack2 = trackWithSecondMaximumRho->getECLCluster();
837 if (eclTrack2) {
838 hasClus++;
840 }
841 const ROOT::Math::PxPyPzEVector V4p1 = PCmsLabTransform::labToCms(trackWithMaximumRho->get4Vector());
842 const ROOT::Math::PxPyPzEVector V4p2 = PCmsLabTransform::labToCms(trackWithSecondMaximumRho->get4Vector());
843
844 const ROOT::Math::PxPyPzEVector V4pSum = V4p1 + V4p2;
845 const double mSum = V4pSum.M();
846
847 const double thetaSumCMS = (V4p1.Theta() + V4p2.Theta()) * TMath::RadToDeg();
848 const double phi1CMS = V4p1.Phi() * TMath::RadToDeg();
849 const double phi2CMS = V4p2.Phi() * TMath::RadToDeg();
850
851 double diffPhi = phi1CMS - phi2CMS;
852 if (std::abs(diffPhi) > 180) {
853 if (diffPhi > 180) {
854 diffPhi = diffPhi - 2 * 180;
855 } else {
856 diffPhi = 2 * 180 + diffPhi;
857 }
858 }
859 const double delThetaCMS = std::abs(std::abs(thetaSumCMS) - 180);
860 const double delPhiCMS = std::abs(180 - std::abs(diffPhi));
861
862 const bool mumuHighMassCand = chSum == 0 && (mSum > 8. && mSum < 12.) && hasClus > 0 && eclE1 <= 1
863 && eclE2 <= 1 && delThetaCMS < 10 && delPhiCMS < 10;
864
865 if (mumuHighMassCand) mumuHighMass = 1;
866
867 }
868
869 calculationResult["MumuHighM"] = mumuHighMass;
870
871 // BtoCharm skims
872 calculationResult["Bp"] = 0;
873 calculationResult["Bz"] = 0;
874
875 if (m_BpParticles.isValid() && (ntrk_bha >= 3 && Bhabha2Trk == 0)) {
876 calculationResult["Bp"] = m_BpParticles->getListSize() >= 1;
877 }
878
879 if (m_BzParticles.isValid() && (ntrk_bha >= 3 && Bhabha2Trk == 0)) {
880 calculationResult["Bz"] = m_BzParticles->getListSize() >= 1;
881 }
882
883}
double getDedxNoSat() const
Get dE/dx truncated mean without the saturation correction for this track.
Class to provide momentum-related information from ECLClusters.
const ROOT::Math::PxPyPzEVector Get4MomentumFromCluster(const ECLCluster *cluster, ECLCluster::EHypothesisBit hypo)
Returns four momentum vector.
static const ChargedStable pion
charged pion particle
Definition Const.h:661
ECL cluster data.
Definition ECLCluster.h:27
bool hasHypothesis(EHypothesisBit bitmask) const
Return if specific hypothesis bit is set.
Definition ECLCluster.h:351
double getEnergy(EHypothesisBit hypothesis) const
Return Energy (GeV).
Definition ECLCluster.cc:23
@ c_nPhotons
CR is split into n photons (N1)
Definition ECLCluster.h:41
unsigned short getNHits() const
Get the total Number of CDC hits in the fit.
ROOT::Math::PxPyPzEVector getBeamFourMomentum() const
Returns LAB four-momentum of e+e-, i.e.
static ROOT::Math::PxPyPzMVector labToCms(const ROOT::Math::PxPyPzMVector &vec)
Transforms Lorentz vector into CM System.
bool isAvailable(Const::PIDDetectorSet set=Const::PIDDetectorSet::set()) const
Check whether PID information is available for at least one of the detectors in a given set.
Class to store reconstructed particles.
Definition Particle.h:76
const Track * getTrack() const
Returns the pointer to the Track object that was used to create this Particle (ParticleType == c_Trac...
Definition Particle.cc:916
const ECLCluster * getECLCluster() const
Returns the pointer to the ECLCluster object that was used to create this Particle (if ParticleType =...
Definition Particle.cc:976
const PIDLikelihood * getPIDLikelihood() const
Returns the pointer to the PIDLikelihood object that is related to the Track, which was used to creat...
Definition Particle.cc:947
double getCharge(void) const
Returns particle charge.
Definition Particle.cc:653
ROOT::Math::PxPyPzEVector get4Vector() const
Returns Lorentz vector.
Definition Particle.h:567
ROOT::Math::XYZVector getMomentum() const
Returns momentum vector.
Definition Particle.h:580
double getMomentumMagnitude() const
Returns momentum magnitude.
Definition Particle.h:589
double getP() const
Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)
Definition Particle.h:598
const Particle * getDaughter(unsigned i) const
Returns a pointer to the i-th daughter particle.
Definition Particle.cc:662
double getExtraInfo(const std::string &name) const
Return given value if set.
Definition Particle.cc:1374
double getMass() const
Returns invariant mass (= nominal for FS particles)
Definition Particle.h:527
static const ReferenceFrame & GetCurrent()
Get current rest frame.
StoreObjPtr< ParticleList > m_pionParticles
Internal storage of the tracks as particles.
StoreObjPtr< ParticleList > m_gammaParticles
Internal storage of the ECL clusters as particles.
StoreObjPtr< ParticleList > m_BzParticles
Internal storage of the B0's.
StoreObjPtr< ParticleList > m_LambdaParticles
Internal storage of the Lambda0's.
StoreObjPtr< ParticleList > m_pionHadParticles
Internal storage of the tracks as particles (definition for hadronb).
void requireStoreArrays() override
Require the particle list. We do not need more here.
void doCalculation(SoftwareTriggerObject &calculationResult) override
Actually write out the variables into the map.
SkimSampleCalculator()
Set the default names for the store object particle lists.
StoreObjPtr< ParticleList > m_KsParticles
Internal storage of the K_S0's.
std::string m_filterL1TrgNN
HLT filter line for the TRG skim.
StoreObjPtr< ParticleList > m_DstParticles
Internal storage of the D*'s.
StoreObjPtr< ParticleList > m_pionTauParticles
Internal storage of the tracks as particles (definition for tau skims).
StoreObjPtr< ParticleList > m_BpParticles
Internal storage of the B+'s.
StoreObjPtr< ParticleList > m_offIpParticles
Internal storage of the tracks for alignment calibration.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
bool isValid() const
Check whether the array was registered.
Definition StoreArray.h:288
int getEntries() const
Get the number of objects in the array.
Definition StoreArray.h:216
iterator end()
Return iterator to last entry +1.
Definition StoreArray.h:320
iterator begin()
Return iterator to first entry.
Definition StoreArray.h:318
bool isValid() const
Check whether the object was created.
ROOT::Math::PxPyPzEVector get4Momentum() const
Getter for the 4Momentum at the closest approach of the track in the r/phi projection.
HitPatternCDC getHitPatternCDC() const
Getter for the hit pattern in the CDC;.
Abstract base class for different kinds of events.