Belle II Software development
KLMClusterVariables.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/* Own header. */
10#include <analysis/variables/KLMClusterVariables.h>
11
12// include VariableManager
13#include <analysis/VariableManager/Manager.h>
14
15/* Analysis headers. */
16#include <analysis/dataobjects/Particle.h>
17#include <analysis/utility/PCmsLabTransform.h>
18
19/* Basf2 headers. */
20#include <framework/datastore/StoreArray.h>
21#include <mdst/dataobjects/ECLCluster.h>
22#include <mdst/dataobjects/KlId.h>
23#include <mdst/dataobjects/KLMCluster.h>
24#include <mdst/dataobjects/Track.h>
25#include <mdst/dataobjects/TrackFitResult.h>
26
27using namespace std;
28
29namespace Belle2::Variable {
30
31 double klmClusterKlId(const Particle* particle)
32 {
33 const KLMCluster* cluster = particle->getKLMCluster();
34 if (!cluster) {
35 return Const::doubleNaN;
36 }
37 const KlId* klid = cluster->getRelatedTo<KlId>();
38 if (!klid) {
39 return Const::doubleNaN;
40 }
41 return klid->getKlId();
42 }
43
44 int klmClusterBelleTrackFlag(const Particle* particle)
45 {
46 const float angle = 0.24;
47 const KLMCluster* cluster = particle->getKLMCluster();
48 if (!cluster) {
49 return 0;
50 }
51 const B2Vector3D& pos = cluster->getClusterPosition();
52 StoreArray<TrackFitResult> tracks;
53 for (const TrackFitResult& track : tracks) {
54 const B2Vector3D& trackPos = track.getPosition();
55 if (trackPos.Angle(pos) < angle) {
56 return 1;
57 }
58 }
59 return 0;
60 }
61
62 int klmClusterBelleECLFlag(const Particle* particle)
63 {
64 const float angle = 0.24;
65 const KLMCluster* klmCluster = particle->getKLMCluster();
66 if (!klmCluster) {
67 return 0;
68 }
69 const B2Vector3D& klmClusterPos = klmCluster->getClusterPosition();
70 StoreArray<ECLCluster> eclClusters;
71 for (const ECLCluster& eclCluster : eclClusters) {
72 const B2Vector3D& eclClusterPos = eclCluster.getClusterPosition();
73 if (eclClusterPos.Angle(klmClusterPos) < angle) {
74 return 1;
75 }
76 }
77 return 0;
78 }
79
80 double klmClusterTiming(const Particle* particle)
81 {
82 const KLMCluster* cluster = particle->getKLMCluster();
83 if (!cluster) {
84 return Const::doubleNaN;
85 }
86 return cluster->getTime();
87 }
88
89
90 double klmClusterPositionX(const Particle* particle)
91 {
92 const KLMCluster* cluster = particle->getKLMCluster();
93 if (!cluster) {
94 return Const::doubleNaN;
95 }
96 return cluster->getClusterPosition().X();
97 }
98
99
100 double klmClusterPositionY(const Particle* particle)
101 {
102 const KLMCluster* cluster = particle->getKLMCluster();
103 if (!cluster) {
104 return Const::doubleNaN;
105 }
106 return cluster->getClusterPosition().Y();
107 }
108
109
110 double klmClusterPositionZ(const Particle* particle)
111 {
112 const KLMCluster* cluster = particle->getKLMCluster();
113 if (!cluster) {
114 return Const::doubleNaN;
115 }
116 return cluster->getClusterPosition().Z();
117 }
118
119
120 double klmClusterInnermostLayer(const Particle* particle)
121 {
122 const KLMCluster* cluster = particle->getKLMCluster();
123 if (!cluster) {
124 return Const::doubleNaN;
125 }
126 return cluster->getInnermostLayer();
127 }
128
129
130 double klmClusterLayers(const Particle* particle)
131 {
132 const KLMCluster* cluster = particle->getKLMCluster();
133 if (!cluster) {
134 return Const::doubleNaN;
135 }
136 return cluster->getLayers();
137 }
138
139 double klmClusterEnergy(const Particle* particle)
140 {
141 const KLMCluster* cluster = particle->getKLMCluster();
142 if (!cluster) {
143 return Const::doubleNaN;
144 }
145 return cluster->getEnergy();
146 }
147
148 double klmClusterMomentum(const Particle* particle)
149 {
150 const KLMCluster* cluster = particle->getKLMCluster();
151 if (!cluster) {
152 return Const::doubleNaN;
153 }
154 return cluster->getMomentumMag();
155 }
156
157 double klmClusterIsBKLM(const Particle* particle)
158 {
159 const KLMCluster* cluster = particle->getKLMCluster();
160 if (!cluster) {
161 return Const::doubleNaN;
162 }
163 float clusterZ = cluster->getClusterPosition().Z();
164 if ((clusterZ > -180) && (clusterZ < 275)) {
165 return 1;
166 }
167 return 0;
168 }
169
170 double klmClusterIsEKLM(const Particle* particle)
171 {
172 const KLMCluster* cluster = particle->getKLMCluster();
173 if (!cluster) {
174 return Const::doubleNaN;
175 }
176 float clusterZ = cluster->getClusterPosition().Z();
177 if ((clusterZ < -180) || (clusterZ > 275)) {
178 return 1;
179 }
180 return 0;
181 }
182
183 double klmClusterIsForwardEKLM(const Particle* particle)
184 {
185 const KLMCluster* cluster = particle->getKLMCluster();
186 if (!cluster) {
187 return Const::doubleNaN;
188 }
189 float clusterZ = cluster->getClusterPosition().Z();
190 if (clusterZ > 275) {
191 return 1;
192 }
193 return 0;
194 }
195
196 double klmClusterIsBackwardEKLM(const Particle* particle)
197 {
198 const KLMCluster* cluster = particle->getKLMCluster();
199 if (!cluster) {
200 return Const::doubleNaN;
201 }
202 float clusterZ = cluster->getClusterPosition().Z();
203 if (clusterZ < -180) {
204 return 1;
205 }
206 return 0;
207 }
208
209 double klmClusterTheta(const Particle* particle)
210 {
211 const KLMCluster* cluster = particle->getKLMCluster();
212 if (!cluster) {
213 return Const::doubleNaN;
214 }
215 return cluster->getClusterPosition().Theta();
216 }
217
218 double klmClusterPhi(const Particle* particle)
219 {
220 const KLMCluster* cluster = particle->getKLMCluster();
221 if (!cluster) {
222 return Const::doubleNaN;
223 }
224 return cluster->getClusterPosition().Phi();
225 }
226
227 double maximumKLMAngleCMS(const Particle* particle)
228 {
229 // check there actually are KLM clusters in the event
230 StoreArray<KLMCluster> clusters;
231 if (clusters.getEntries() == 0) return Const::doubleNaN;
232
233 // get the input particle's vector momentum in the CMS frame
234 PCmsLabTransform T;
235 const B2Vector3D pCms = (T.rotateLabToCms() * particle->get4Vector()).Vect();
236
237 // find the KLM cluster with the largest angle
238 double maxAngle = 0.0;
239 for (int iKLM = 0; iKLM < clusters.getEntries(); iKLM++) {
240 const B2Vector3D clusterMomentumCms = (T.rotateLabToCms() * clusters[iKLM]->getMomentum()).Vect();
241 double angle = pCms.Angle(clusterMomentumCms);
242 if (angle > maxAngle) maxAngle = angle;
243 }
244 return maxAngle;
245 }
246
247 double nKLMClusterTrackMatches(const Particle* particle)
248 {
249 const KLMCluster* cluster = particle->getKLMCluster();
250 if (!cluster)
251 return Const::doubleNaN;
252 size_t out = cluster->getRelationsFrom<Track>().size();
253 return double(out);
254 }
255
256 double nMatchedKLMClusters(const Particle* particle)
257 {
258 Belle2::Particle::EParticleSourceObject particleSource = particle->getParticleSource();
259 if (particleSource == Particle::EParticleSourceObject::c_Track) {
260 return particle->getTrack()->getRelationsTo<KLMCluster>().size();
261 } else if (particleSource == Particle::EParticleSourceObject::c_ECLCluster) {
262 return particle->getECLCluster()->getRelationsTo<KLMCluster>().size();
263 } else {
264 return Const::doubleNaN;
265 }
266 }
267
268 double nKLMClusterECLClusterMatches(const Particle* particle)
269 {
270 const KLMCluster* cluster = particle->getKLMCluster();
271 if (!cluster)
272 return Const::doubleNaN;
273 size_t out = cluster->getRelationsFrom<ECLCluster>().size();
274 return double(out);
275 }
276
277 double klmClusterTrackDistance(const Particle* particle)
278 {
279 const KLMCluster* cluster = particle->getKLMCluster();
280 if (!cluster)
281 return Const::doubleNaN;
282 auto trackWithWeight = cluster->getRelatedFromWithWeight<Track>();
283 if (!trackWithWeight.first)
284 return Const::doubleNaN;
285 return 1. / trackWithWeight.second;
286 }
287
288 VARIABLE_GROUP("KLM Cluster and KlongID");
289
290 REGISTER_VARIABLE("klmClusterKlId", klmClusterKlId,
291 "Returns the KlId classifier output associated to the KLMCluster.");
292 REGISTER_VARIABLE("klmClusterBelleTrackFlag", klmClusterBelleTrackFlag,
293 "Returns the Belle-style Track flag.");
294 REGISTER_VARIABLE("klmClusterBelleECLFlag", klmClusterBelleECLFlag,
295 "Returns the Belle-style ECL flag.");
296 REGISTER_VARIABLE("klmClusterTiming", klmClusterTiming, R"DOC(
297Returns the timing information of the associated KLMCluster.
298
299)DOC","ns");
300 REGISTER_VARIABLE("klmClusterPositionX", klmClusterPositionX, R"DOC(
301Returns the :math:`x` position of the associated KLMCluster.
302
303)DOC","cm");
304 REGISTER_VARIABLE("klmClusterPositionY", klmClusterPositionY, R"DOC(
305Returns the :math:`y` position of the associated KLMCluster.
306
307)DOC","cm");
308 REGISTER_VARIABLE("klmClusterPositionZ", klmClusterPositionZ, R"DOC(
309Returns the :math:`z` position of the associated KLMCluster.
310
311)DOC","cm");
312 REGISTER_VARIABLE("klmClusterInnermostLayer", klmClusterInnermostLayer,
313 "Returns the number of the innermost KLM layer with a 2-dimensional hit of the associated KLMCluster.");
314 REGISTER_VARIABLE("klmClusterLayers", klmClusterLayers,
315 "Returns the number of KLM layers with 2-dimensional hits of the associated KLMCluster.");
316 REGISTER_VARIABLE("klmClusterEnergy", klmClusterEnergy, R"DOC(
317Returns the energy of the associated KLMCluster.
318
319.. warning::
320 This variable returns an approximation of the energy: it uses :b2:var:`klmClusterMomentum` as momentum and the hypothesis that the KLMCluster is originated by a :math:`K_{L}^0`
321 (:math:`E_{\text{KLM}} = \sqrt{M_{K^0_L}^2 + p_{\text{KLM}}^2}`, where :math:`E_{\text{KLM}}` is this variable, :math:`M_{K^0_L}` is the :math:`K^0_L` mass and :math:`p_{\text{KLM}}` is :b2:var:`klmClusterMomentum`).
322 It should be used with caution, and may not be physically meaningful, especially for :math:`n` candidates.
323
324)DOC","GeV");
325 REGISTER_VARIABLE("klmClusterMomentum", klmClusterMomentum, R"DOC(
326Returns the momentum magnitude of the associated KLMCluster.
327
328.. warning::
329 This variable returns an approximation of the momentum, since it is proportional to :b2:var:`klmClusterLayers`
330 (:math:`p_{\text{KLM}} = 0.215 \cdot N_{\text{layers}}`, where :math:`p_{\text{KLM}}` is this variable and :math:`N_{\text{layers}}` is :b2:var:`klmClusterLayers`).
331 It should be used with caution, and may not be physically meaningful.
332
333)DOC","GeV/c");
334 REGISTER_VARIABLE("klmClusterIsBKLM", klmClusterIsBKLM,
335 "Returns 1 if the associated KLMCluster is in barrel KLM.");
336 REGISTER_VARIABLE("klmClusterIsEKLM", klmClusterIsEKLM,
337 "Returns 1 if the associated KLMCluster is in endcap KLM.");
338 REGISTER_VARIABLE("klmClusterIsForwardEKLM", klmClusterIsForwardEKLM,
339 "Returns 1 if the associated KLMCluster is in forward endcap KLM.");
340 REGISTER_VARIABLE("klmClusterIsBackwardEKLM", klmClusterIsBackwardEKLM,
341 "Returns 1 if the associated KLMCluster is in backward endcap KLM.");
342 REGISTER_VARIABLE("klmClusterTheta", klmClusterTheta, R"DOC(
343Returns the polar (:math:`\theta`) angle of the associated KLMCluster.
344
345)DOC","rad");
346 REGISTER_VARIABLE("klmClusterPhi", klmClusterPhi, R"DOC(
347Returns the azimuthal (:math:`\phi`) angle of the associated KLMCluster.
348
349)DOC","rad");
350 REGISTER_VARIABLE("maximumKLMAngleCMS", maximumKLMAngleCMS ,
351 "Returns the maximum angle in the CMS frame between the Particle and all KLMClusters in the event.\n\n","rad");
352 REGISTER_VARIABLE("nKLMClusterTrackMatches", nKLMClusterTrackMatches, R"DOC(
353Returns the number of Tracks matched to the KLMCluster associated to this Particle. This variable can return a number greater than 0 for :math:`K_{L}^0` or :math:`n` candidates originating from KLMClusters and returns NaN for Particles with no KLMClusters associated.
354)DOC");
355 REGISTER_VARIABLE("nMatchedKLMClusters", nMatchedKLMClusters, R"DOC(
356 Returns the number of KLMClusters matched to the particle. It only works for
357 Particles created either from Tracks or from ECLCluster, while it returns NaN
358 for :math:`K_{L}^0` or :math:`n` candidates originating from KLMClusters.
359 )DOC");
360 REGISTER_VARIABLE("nKLMClusterECLClusterMatches", nKLMClusterECLClusterMatches, R"DOC(
361 Returns the number of ECLClusters matched to the KLMCluster associated to this Particle.
362 )DOC");
363 REGISTER_VARIABLE("klmClusterTrackDistance", klmClusterTrackDistance,
364 "Returns the distance between the Track and the KLMCluster associated to this Particle. This variable returns NaN if there is no Track-to-KLMCluster relationship.\n\n",
365 "cm");
366
367}
static const double doubleNaN
quiet_NaN
Definition: Const.h:703
EParticleSourceObject
particle source enumerators
Definition: Particle.h:82
B2Vector3< double > B2Vector3D
typedef for common usage with double
Definition: B2Vector3.h:516
STL namespace.