Belle II Software development
DataWriterModule.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
10
11#include <reconstruction/modules/KlId/DataWriter/DataWriterModule.h>
12#include <mdst/dataobjects/KlId.h>
13
14#include <framework/datastore/StoreArray.h>
15
16#include <tracking/dataobjects/TrackClusterSeparation.h>
17
18#include <analysis/ClusterUtility/ClusterUtils.h>
19
20#include <mdst/dataobjects/MCParticle.h>
21#include <mdst/dataobjects/KLMCluster.h>
22#include <mdst/dataobjects/ECLCluster.h>
23
24#include <TTree.h>
25#include <TFile.h>
26#include <utility>
27#include <Math/VectorUtil.h>
28
29#include <reconstruction/modules/KlId/KLMExpert/KlId.h>
30
31using namespace std;
32using namespace Belle2;
33using namespace Belle2::KlongId;
34
35
36// --------------------------------------Module----------------------------------------------
37
38
39REG_MODULE(DataWriter);
40
42{
43 setDescription("Used to write flat ntuple for KlId classifier trainings for both ECL and KLM KlID. Output is a root file.");
44
45 addParam("outPath", m_outPath, "Output path - where you want your root files to be placed.", m_outPath);
46 addParam("useKLM", m_useKLM, "Write KLM data.", m_useKLM);
47 addParam("useECL", m_useECL, "Write ECL data.", m_useECL);
48}
49
50
51
55
56
58{
59 // require existence of necessary datastore obj
60
61 m_eclClusters.isRequired();
62 m_klmClusters.isRequired();
63 m_mcParticles.isRequired();
64
65 m_klmClusters.requireRelationTo(m_mcParticles);
66 m_klmClusters.registerRelationTo(m_eclClusters);
67
68 m_f = new TFile(m_outPath.c_str(), "recreate");
69
70 m_treeECLhadron = new TTree("ECLdataHadron", "ECLdataHadron");
71 m_treeECLgamma = new TTree("ECLdataGamma", "ECLdataGamma");
72
73 // KLM
74 if (m_useKLM) {
75 m_treeKLM = new TTree("KLMdata", "KLMdata");
76 m_treeKLM -> Branch("KLMMCMom", & m_KLMMCMom);
77 m_treeKLM -> Branch("KLMMCPhi", & m_KLMMCPhi);
78 m_treeKLM -> Branch("KLMMCTheta", & m_KLMMCTheta);
79 m_treeKLM -> Branch("KLMMom", & m_KLMMom);
80 m_treeKLM -> Branch("KLMPhi", & m_KLMPhi);
81 m_treeKLM -> Branch("KLMTheta", & m_KLMTheta);
82 m_treeKLM -> Branch("KLMMCLifetime", & m_KLMMCLifetime);
83 m_treeKLM -> Branch("KLMMCPDG", & m_KLMMCPDG);
84 m_treeKLM -> Branch("KLMMCPrimaryPDG", & m_KLMMCPrimaryPDG);
85 m_treeKLM -> Branch("KLMMCStatus", & m_KLMMCStatus);
86 m_treeKLM -> Branch("KLMnCluster", & m_KLMnCluster);
87 m_treeKLM -> Branch("KLMnLayer", & m_KLMnLayer);
88 m_treeKLM -> Branch("KLMnInnermostlayer", & m_KLMnInnermostLayer);
89 m_treeKLM -> Branch("KLMglobalZ", & m_KLMglobalZ);
90 m_treeKLM -> Branch("KLMtime", & m_KLMtime);
91 m_treeKLM -> Branch("KLMinvM", & m_KLMinvM);
92 m_treeKLM -> Branch("KLMTruth", & m_KLMTruth);
93 m_treeKLM -> Branch("KLMdistToNextCl", & m_KLMnextCluster);
94 m_treeKLM -> Branch("KLMenergy", & m_KLMenergy);
95 m_treeKLM -> Branch("KLMaverageInterClusterDist", & m_KLMavInterClusterDist);
96 m_treeKLM -> Branch("KLMhitDepth", & m_KLMhitDepth);
97
98 m_treeKLM -> Branch("KLMECLHypo", & m_KLMECLHypo);
99 m_treeKLM -> Branch("KLMECLZMVA", & m_KLMECLZMVA);
100 m_treeKLM -> Branch("KLMECLZ40", & m_KLMECLZ40);
101 m_treeKLM -> Branch("KLMECLZ51", & m_KLMECLZ51);
102 m_treeKLM -> Branch("KLMECLUncertaintyPhi", & m_KLMECLUncertaintyPhi);
103 m_treeKLM -> Branch("KLMECLUncertaintyTheta", & m_KLMECLUncertaintyTheta);
104 m_treeKLM -> Branch("KLMdistToNextECL", & m_KLMECLDist);
105 m_treeKLM -> Branch("KLMtrackToECL", & m_KLMtrackToECL);
106 m_treeKLM -> Branch("KLMECLEerror", & m_KLMECLEerror);
107 m_treeKLM -> Branch("KLMECLenergy", & m_KLMECLE);
108 m_treeKLM -> Branch("KLMECLE9oE25", & m_KLMECLE9oE25);
109 m_treeKLM -> Branch("KLMECLtiming", & m_KLMECLTiming);
110 m_treeKLM -> Branch("KLMECLTerror", & m_KLMECLTerror);
111 m_treeKLM -> Branch("KLMECLdeltaL", & m_KLMECLdeltaL);
112 m_treeKLM -> Branch("KLMECLmintrackDist", & m_KLMECLminTrackDist);
113 m_treeKLM -> Branch("KLMTrackSepDist", & m_KLMTrackSepDist);
114 m_treeKLM -> Branch("KLMTrackSepAngle", & m_KLMTrackSepAngle);
115 m_treeKLM -> Branch("KLMInitialtrackSepAngle", & m_KLMInitialTrackSepAngle);
116 m_treeKLM -> Branch("KLMTrackRotationAngle", & m_KLMTrackRotationAngle);
117 m_treeKLM -> Branch("KLMTrackClusterSepAngle", & m_KLMTrackClusterSepAngle);
118 m_treeKLM -> Branch("isBeamBKG", & m_isBeamBKG);
119 m_treeKLM -> Branch("KLMKlId", & m_KLMKLid);
120 m_treeKLM -> Branch("KLMAngleToMC", & m_KLMAngleToMC);
121 m_treeKLM -> Branch("KLMMCWeight", & m_KLMMCWeight);
122 m_treeKLM -> Branch("KLMtrackFlag", & m_KLMtrackFlag);
123 m_treeKLM -> Branch("KLMeclFlag", & m_KLMeclFlag);
124 m_treeKLM -> Branch("isSignal", & m_isSignal);
125 }//useKLM
126
127 //ECL
128 if (m_useECL) {
129 m_treeECLhadron = new TTree("ECLdataHadron", "ECLdataHadron");
130 m_treeECLgamma = new TTree("ECLdataGamma", "ECLdataGamma");
131
132 m_treeECLhadron -> Branch("ECLMCMom", & m_ECLMCMom);
133 m_treeECLhadron -> Branch("ECLMCPhi", & m_ECLMCPhi);
134 m_treeECLhadron -> Branch("ECLMCLifetime", & m_ECLMCLifetime);
135 m_treeECLhadron -> Branch("ECLMCPDG", & m_ECLMCPDG);
136 m_treeECLhadron -> Branch("ECLMCTheta", & m_ECLMCTheta);
137 m_treeECLhadron -> Branch("ECLMCLifetime", & m_ECLMCLifetime);
138 m_treeECLhadron -> Branch("ECLMCPDG", & m_ECLMCPDG);
139 m_treeECLhadron -> Branch("ECLMCStatus", & m_ECLMCStatus);
140 m_treeECLhadron -> Branch("ECLMCPrimaryPDG", & m_ECLMCPrimaryPDG);
141 m_treeECLhadron -> Branch("ECLUncertaintyEnergy", & m_ECLUncertaintyEnergy);
142 m_treeECLhadron -> Branch("ECLUncertaintyTheta", & m_ECLUncertaintyTheta);
143 m_treeECLhadron -> Branch("ECLUncertaintyPhi", & m_ECLUncertaintyPhi);
144 m_treeECLhadron -> Branch("ECLMom", & m_ECLMom);
145 m_treeECLhadron -> Branch("ECLPhi", & m_ECLPhi);
146 m_treeECLhadron -> Branch("ECLTheta", & m_ECLTheta);
147 m_treeECLhadron -> Branch("ECLZ", & m_ECLZ);
148 m_treeECLhadron -> Branch("ECLenergy", & m_ECLE);
149 m_treeECLhadron -> Branch("ECLE9oE25", & m_ECLE9oE25);
150 m_treeECLhadron -> Branch("ECLtiming", & m_ECLTiming);
151 m_treeECLhadron -> Branch("ECLR", & m_ECLR);
152 m_treeECLhadron -> Branch("ECLTruth", & m_ECLTruth);
153 m_treeECLhadron -> Branch("ECLZ51", & m_ECLZ51);
154 m_treeECLhadron -> Branch("ECLZ40", & m_ECLZ40);
155 m_treeECLhadron -> Branch("ECLE1oE9", & m_ECLE1oE9);
156 m_treeECLhadron -> Branch("ECL2ndMom", & m_ECL2ndMom);
157 m_treeECLhadron -> Branch("ECLnumChrystals", & m_ECLnumChrystals);
158 m_treeECLhadron -> Branch("ECLLAT", & m_ECLLAT);
159 m_treeECLhadron -> Branch("ECLZMVA", & m_ECLZMVA);
160 m_treeECLhadron -> Branch("ECLKlId", & m_ECLKLid);
161 m_treeECLhadron -> Branch("ECLdeltaL", & m_ECLdeltaL);
162 m_treeECLhadron -> Branch("ECLmintrackDist", & m_ECLminTrkDistance);
163 m_treeECLhadron -> Branch("isBeamBKG", & m_isBeamBKG);
164 m_treeECLhadron -> Branch("ECLMCWeight", & m_ECLMCWeight);
165 m_treeECLhadron -> Branch("isSignal", & m_isSignal);
166
167
168 m_treeECLgamma -> Branch("ECLMCMom", & m_ECLMCMom);
169 m_treeECLgamma -> Branch("ECLMCPhi", & m_ECLMCPhi);
170 m_treeECLgamma -> Branch("ECLMCTheta", & m_ECLMCTheta);
171 m_treeECLgamma -> Branch("ECLMCLifetime", & m_ECLMCLifetime);
172 m_treeECLgamma -> Branch("ECLMCPDG", & m_ECLMCPDG);
173 m_treeECLgamma -> Branch("ECLMCStatus", & m_ECLMCStatus);
174 m_treeECLgamma -> Branch("ECLMCPrimaryPDG", & m_ECLMCPrimaryPDG);
175 m_treeECLgamma -> Branch("ECLUncertaintyEnergy", & m_ECLUncertaintyEnergy);
176 m_treeECLgamma -> Branch("ECLUncertaintyTheta", & m_ECLUncertaintyTheta);
177 m_treeECLgamma -> Branch("ECLUncertaintyPhi", & m_ECLUncertaintyPhi);
178 m_treeECLgamma -> Branch("ECLMom", & m_ECLMom);
179 m_treeECLgamma -> Branch("ECLPhi", & m_ECLPhi);
180 m_treeECLgamma -> Branch("ECLTheta", & m_ECLTheta);
181 m_treeECLgamma -> Branch("ECLZ", & m_ECLZ);
182 m_treeECLgamma -> Branch("ECLenergy", & m_ECLE);
183 m_treeECLgamma -> Branch("ECLE9oE25", & m_ECLE9oE25);
184 m_treeECLgamma -> Branch("ECLtiming", & m_ECLTiming);
185 m_treeECLgamma -> Branch("ECLR", & m_ECLR);
186 m_treeECLgamma -> Branch("ECLTruth", & m_ECLTruth);
187 m_treeECLgamma -> Branch("ECLZ51", & m_ECLZ51);
188 m_treeECLgamma -> Branch("ECLZ40", & m_ECLZ40);
189 m_treeECLgamma -> Branch("ECLE1oE9", & m_ECLE1oE9);
190 m_treeECLgamma -> Branch("ECL2ndMom", & m_ECL2ndMom);
191 m_treeECLgamma -> Branch("ECLnumChrystals", & m_ECLnumChrystals);
192 m_treeECLgamma -> Branch("ECLLAT", & m_ECLLAT);
193 m_treeECLgamma -> Branch("ECLZMVA", & m_ECLZMVA);
194 m_treeECLgamma -> Branch("ECLKlId", & m_ECLKLid);
195 m_treeECLgamma -> Branch("ECLdeltaL", & m_ECLdeltaL);
196 m_treeECLgamma -> Branch("ECLmintrackDist", & m_ECLminTrkDistance);
197 m_treeECLgamma -> Branch("isBeamBKG", & m_isBeamBKG);
198 m_treeECLgamma -> Branch("ECLMCWeight", & m_ECLMCWeight);
199 m_treeECLgamma -> Branch("isSignal", & m_isSignal);
200 }//useECL
201}//init
202
203
207
211
213{
214 // Use the neutralHadron hypothesis for the ECL
216
217 for (const KLMCluster& cluster : m_klmClusters) {
218
219 if (!m_useKLM) {continue;}
220
221 const ROOT::Math::XYZVector& clusterPos = cluster.getClusterPosition();
222
223 m_KLMPhi = clusterPos.Phi();
224 m_KLMTheta = clusterPos.Theta();
225
226 m_KLMglobalZ = clusterPos.Z();
227 m_KLMnCluster = m_klmClusters.getEntries();
228 m_KLMnLayer = cluster.getLayers();
229 m_KLMnInnermostLayer = cluster.getInnermostLayer();
230 m_KLMtime = cluster.getTime();
231 m_KLMinvM = cluster.getMomentum().M2();
232 m_KLMenergy = cluster.getEnergy();
233 m_KLMhitDepth = cluster.getClusterPosition().R();
234 m_KLMtrackFlag = cluster.getAssociatedTrackFlag();
235 m_KLMeclFlag = cluster.getAssociatedEclClusterFlag();
236
237 m_KLMTrackSepDist = -999;
238 m_KLMTrackSepAngle = -999;
242 auto trackSeperations = cluster.getRelationsTo<TrackClusterSeparation>();
243 TrackClusterSeparation* trackSep;
244 float best_dist = 100000000;
245 for (auto trackSeperation : trackSeperations) {
246 float dist = trackSeperation.getDistance();
247 if (dist < best_dist) {
248 best_dist = dist;
249 trackSep = &trackSeperation;
250 m_KLMTrackSepDist = trackSep->getDistance();
251 m_KLMTrackSepAngle = trackSep->getTrackClusterAngle();
252 m_KLMInitialTrackSepAngle = trackSep->getTrackClusterInitialSeparationAngle();
253 m_KLMTrackRotationAngle = trackSep->getTrackRotationAngle();
254 m_KLMTrackClusterSepAngle = trackSep->getTrackClusterSeparationAngle();
255 }
256 }
257
258
259 if (isnan(m_KLMglobalZ)) { m_KLMglobalZ = -999;}
260 if (isnan(m_KLMnCluster)) { m_KLMnCluster = -999;}
261 if (isnan(m_KLMnLayer)) { m_KLMnLayer = -999;}
262 if (isnan(m_KLMnInnermostLayer)) { m_KLMnInnermostLayer = -999;}
263 if (isnan(m_KLMtime)) { m_KLMtime = -999;}
264 if (isnan(m_KLMinvM)) { m_KLMinvM = -999;}
265 if (isnan(m_KLMenergy)) { m_KLMenergy = -999;}
266 if (isnan(m_KLMhitDepth)) { m_KLMhitDepth = -999;}
267 if (isnan(m_KLMTrackSepDist)) { m_KLMTrackSepDist = -999;}
268 if (isnan(m_KLMTrackSepAngle)) { m_KLMTrackSepAngle = -999;}
272
273 // findClosestECLCluster with the c_neutralHadron hypothesis
274 pair<ECLCluster*, double> closestECLAndDist = findClosestECLCluster(clusterPos, eclHypothesis);
275 ECLCluster* closestECLCluster = get<0>(closestECLAndDist);
276 m_KLMECLDist = get<1>(closestECLAndDist);
277
278 if (!(closestECLCluster == nullptr)) {
279 m_KLMECLE = closestECLCluster->getEnergy(eclHypothesis);
280 m_KLMECLE9oE25 = closestECLCluster->getE9oE21();
281 m_KLMECLEerror = closestECLCluster->getUncertaintyEnergy();
282 m_KLMECLTerror = closestECLCluster->getDeltaTime99();
283 m_KLMECLdeltaL = closestECLCluster->getDeltaL();;
284 m_KLMECLminTrackDist = closestECLCluster->getMinTrkDistance();
285 m_KLMECLTiming = closestECLCluster->getTime();
286 m_KLMECLZMVA = closestECLCluster->getZernikeMVA();
287 m_KLMECLZ40 = closestECLCluster->getAbsZernike40();
288 m_KLMECLZ51 = closestECLCluster->getAbsZernike51();
289 m_KLMECLUncertaintyPhi = closestECLCluster->getUncertaintyPhi();
290 m_KLMECLUncertaintyTheta = closestECLCluster->getUncertaintyTheta();
291 m_KLMECLHypo = closestECLCluster->getHypotheses();
292 } else {
293 m_KLMECLdeltaL = -999;
295 m_KLMECLE = -999;
296 m_KLMECLE9oE25 = -999;
297 m_KLMECLTiming = -999;
298 m_KLMECLTerror = -999;
299 m_KLMECLEerror = -999;
300
301 m_KLMECLHypo = -999;
302 m_KLMECLZMVA = -999;
303 m_KLMECLZ40 = -999;
304 m_KLMECLZ51 = -999;
307 }
308
309 tuple<const KLMCluster*, double, double> closestKLMAndDist = findClosestKLMCluster(clusterPos);
310 m_KLMnextCluster = get<1>(closestKLMAndDist);
311 m_KLMavInterClusterDist = get<2>(closestKLMAndDist);
312
313 const auto mcParticleWeightPair = cluster.getRelatedToWithWeight<MCParticle>();
314 MCParticle* part = mcParticleWeightPair.first;
315
318
319 if (part) {
320 m_KLMMCWeight = mcParticleWeightPair.second;
321 m_KLMAngleToMC = ROOT::Math::VectorUtil::Angle(clusterPos, part->getMomentum());
322 m_KLMMCStatus = part->getStatus();
323 m_KLMMCLifetime = part->getLifetime();
324 m_KLMMCPDG = std::abs(part->getPDG());
326 m_KLMMCMom = part->getMomentum().R();
327 m_KLMMCPhi = part->getMomentum().Phi();
328 m_KLMMCTheta = part->getMomentum().Theta();
329 } else {
330 m_KLMMCWeight = -999;
331 m_KLMAngleToMC = -999;
332 m_KLMMCStatus = -999;
333 m_KLMMCLifetime = -999;
334 m_KLMMCPDG = -999;
335 m_KLMMCPrimaryPDG = -999;
336 m_KLMMCMom = -999;
337 m_KLMMCPhi = -999;
338 m_KLMMCTheta = -999;
339 }
340
341 KlId* klid = cluster.getRelatedTo<KlId>();
342 if (klid) {
343 m_KLMKLid = klid->getKlId();
344 } else {
345 m_KLMKLid = -999;
346 }
347 m_isSignal = isKLMClusterSignal(cluster, 0);
348
349
350 m_treeKLM -> Fill();
351 }// for klmcluster in klmclusters
352
353// --------------- ECL CLUSTERS
354 for (const ECLCluster& cluster : m_eclClusters) {
355
356 if (!m_useECL or !cluster.hasHypothesis(eclHypothesis)) {continue;}
357
358 m_ECLminTrkDistance = cluster.getMinTrkDistance();
359 m_ECLdeltaL = cluster.getDeltaL();
360 m_ECLE = cluster.getEnergy(eclHypothesis);
361 m_ECLE9oE25 = cluster.getE9oE21();
362 m_ECLTiming = cluster.getTime();
363 m_ECLR = cluster.getR();
364 m_ECLEerror = cluster.getUncertaintyEnergy();
365 m_ECLZ51 = cluster.getAbsZernike51();
366 m_ECLZ40 = cluster.getAbsZernike40();
367 m_ECLE1oE9 = cluster.getE1oE9();
368 m_ECL2ndMom = cluster.getSecondMoment();
369 m_ECLnumChrystals = cluster.getNumberOfCrystals();
370 m_ECLLAT = cluster.getLAT();
371 m_ECLZMVA = cluster.getZernikeMVA();
372
373 if (isnan(m_ECLminTrkDistance)) { m_ECLminTrkDistance = -999;}
374 if (isnan(m_ECLdeltaL)) { m_ECLdeltaL = -999;}
375 if (isnan(m_ECLE)) { m_ECLE = -999;}
376 if (isnan(m_ECLE9oE25)) { m_ECLE9oE25 = -999;}
377 if (isnan(m_ECLTiming)) { m_ECLTiming = -999;}
378 if (isnan(m_ECLR)) { m_ECLR = -999;}
379 if (isnan(m_ECLEerror)) { m_ECLEerror = -999;}
380 if (isnan(m_ECLZ40)) { m_ECLZ40 = -999;}
381 if (isnan(m_ECLZ51)) { m_ECLZ51 = -999;}
382 if (isnan(m_ECLE1oE9)) { m_ECLE1oE9 = -999;}
383 if (isnan(m_ECL2ndMom)) { m_ECL2ndMom = -999;}
384 if (isnan(m_ECLnumChrystals)) { m_ECLnumChrystals = -999;}
385 if (isnan(m_ECLLAT)) { m_ECLLAT = -999;}
386 if (isnan(m_ECLZMVA)) { m_ECLZMVA = -999;}
387
388 KlId* klid = cluster.getRelatedTo<KlId>();
389 if (klid) {
390 m_ECLKLid = klid->getKlId();
391 } else {
392 m_ECLKLid = -999;
393 }
394
395 const ROOT::Math::XYZVector& clusterPos = cluster.getClusterPosition();
396
397 m_ECLPhi = clusterPos.Phi();
398 m_ECLTheta = clusterPos.Theta();
399 m_ECLZ = clusterPos.Z();
400
401 ClusterUtils C;
402 m_ECLMom = C.Get4MomentumFromCluster(&cluster, eclHypothesis).Vect().Mag2();
403 m_ECLDeltaTime = cluster.getDeltaTime99();
404
405 m_ECLUncertaintyEnergy = cluster.getUncertaintyEnergy();
406 m_ECLUncertaintyTheta = cluster.getUncertaintyTheta();
407 m_ECLUncertaintyPhi = cluster.getUncertaintyPhi();
408
409// MCParticle* part = cluster.getRelatedTo<MCParticle>();
410 const auto mcParticleWeightPair = cluster.getRelatedToWithWeight<MCParticle>();
411 MCParticle* part = mcParticleWeightPair.first;
412
415 if (part) {
416 m_ECLMCWeight = mcParticleWeightPair.second;
417 m_ECLMCStatus = part->getStatus();
418 m_ECLMCLifetime = part->getLifetime();
419 m_ECLMCPDG = std::abs(part->getPDG());
420
422
423 m_ECLMCMom = part->getMomentum().Mag2();
424 m_ECLMCPhi = part->getMomentum().Phi();
425 m_ECLMCTheta = part->getMomentum().Theta();
426 } else {
427 m_ECLMCWeight = -999;
428 m_ECLMCStatus = -999;
429 m_ECLMCLifetime = -999;
430 m_ECLMCPrimaryPDG = -999;
431 m_ECLMCPDG = -999;
432 m_ECLMCMom = -999;
433 m_ECLMCPhi = -999;
434 m_ECLMCTheta = -999;
435 }
436
438
439 if (cluster.hasHypothesis(ECLCluster::EHypothesisBit::c_neutralHadron)) { m_treeECLhadron -> Fill();}
440 if (cluster.hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)) { m_treeECLgamma -> Fill();}
441 }// for ecl cluster in clusters
442} // event
443
444
446{
447 // close root files
448 m_f -> cd();
449 if (m_useKLM) { m_treeKLM -> Write();}
450 if (m_useECL) {
451 m_treeECLhadron -> Write();
452 m_treeECLgamma -> Write();
453 }
454 m_f -> Close();
455}
456
Class to provide momentum-related information from ECLClusters.
const ROOT::Math::PxPyPzEVector Get4MomentumFromCluster(const ECLCluster *cluster, ECLCluster::EHypothesisBit hypo)
Returns four momentum vector.
Float_t m_KLMnLayer
number of layers hit in KLM cluster
bool m_useKLM
write out KLM data
Float_t m_ECLnumChrystals
number of crystals in the cluster
Float_t m_KLMInitialTrackSepAngle
angular distance from track to cluster at track starting point
Float_t m_KLMTrackSepDist
distance from track separation object
Float_t m_KLMnCluster
varibales to write out.
Float_t m_KLMMCPhi
phi of matched mc particle
Float_t m_ECLMCLifetime
MC particles lifetime.
Float_t m_ECLE1oE9
central crystal devided by 3x3 area with it in its center
Float_t m_KLMTheta
measured theta
Float_t m_KLMECLdeltaL
distance between track entry point and cluster center, might be removed
Float_t m_KLMTrackSepAngle
angular distance from track separation object
Float_t m_ECLMCPrimaryPDG
pdg code of higher order MC particle, a cluster related to a photon that originates from a pi0 decay ...
StoreArray< KLMCluster > m_klmClusters
Store array.
Float_t m_KLMECLminTrackDist
track distance between associated ECL cluster and track extrapolated into ECL
Float_t m_ECLE9oE25
energy of 9/25 chrystall rings (E dispersion shape)
virtual void initialize() override
init
Float_t m_KLMMCLifetime
MC partilces life time.
std::string m_outPath
Output path variable.
Float_t m_ECLLAT
lateral shower shape
virtual void event() override
process event
Float_t m_KLMMom
measured momentum
Float_t m_KLMECLZ51
zernike moment 5,1 of closest ECL cluster
Float_t m_ECLdeltaL
distance between track entrace into cluster and cluster center
Float_t m_KLMMCPDG
pdg code of matched MCparticle
Float_t m_ECLminTrkDistance
more sophisticated distaqnce to track in ECL
Float_t m_ECLMCPDG
pdg code of the MCparticle directly related to the cluster
Float_t m_ECLMCWeight
mc weight
Float_t m_ECLKLid
classifier output
Float_t m_KLMTrackRotationAngle
angle between track at poca and trackbeginning
Float_t m_ECLEerror
uncertainty on E measurement in ECL
Float_t m_ECLR
distance of cluster to IP
Float_t m_ECLMom
measured momentum
Float_t m_ECLMCPhi
MC particle phi; -999 if not MCparticle.
virtual void endRun() override
end run
Float_t m_KLMTruth
target variable for KLM classification
Float_t m_KLMAngleToMC
angle between KLMcluster and Mcparticle
Float_t m_KLMECLE9oE25
E in surrounding 9 crystals divided by surrounding 25 crydtalls.
Float_t m_ECLMCStatus
mc status, seen in detector etc.
Float_t m_ECLE
measured energy
Float_t m_KLMECLZMVA
zernike mva output for closest ECL cluster (based on around 10 z-moments)
virtual void terminate() override
terminate
TTree * m_treeECLhadron
tree containing ntuples for ECL cluster with N2 (hadron hypothesis)
Float_t m_KLMtrackToECL
primitive distance cluster <-> track for associated ECL cluster
Float_t m_ECLZ40
Zernike moment 4,0 see Belle2 note on that.
Float_t m_ECLZ
measured Z-coordinate
Float_t m_ECLUncertaintyPhi
measured uncertainty of phi
Float_t m_ECL2ndMom
second moment, shower shape
Float_t m_KLMECLUncertaintyTheta
theta uncertainty of closest ECL cluster
Float_t m_ECLZMVA
output of a BDT that was fitted on some Zernike Moments on a connected region
TTree * m_treeECLgamma
tree containing ntuples for ECL cluster with N1 (photon hypothesis)
bool m_useECL
write out KLM data
Float_t m_KLMMCStatus
MC particles status.
virtual void beginRun() override
beginn run
Float_t m_isBeamBKG
is beam bkg
Float_t m_KLMMCTheta
theta of matched mc particle
Float_t m_ECLTiming
timing of ECL
Float_t m_KLMMCWeight
mc weight
Float_t m_KLMECLTiming
timing of associated ECL cluster
Float_t m_KLMECLDist
distance associated ECL <-> KLM cluster
Float_t m_KLMglobalZ
global Z position in KLM
StoreArray< ECLCluster > m_eclClusters
Store array.
Float_t m_isSignal
isSignal for the classifier
Float_t m_KLMnInnermostLayer
number of innermost layers hit
Float_t m_ECLTruth
ECL trarget variable.
Float_t m_KLMtime
timing of KLM Cluster
Float_t m_KLMECLTerror
uncertainty on time in associated ECL cluster
Float_t m_KLMinvM
invariant mass calculated from root vector
Float_t m_KLMnextCluster
distance to next KLM cluster
Float_t m_ECLUncertaintyEnergy
measured energy uncertainty
Float_t m_ECLTheta
measured theta
Float_t m_ECLUncertaintyTheta
measured uncertainty on theta
Float_t m_KLMECLHypo
hypotheis id of closest ecl cluster 5: gamma, 6:hadron
StoreArray< MCParticle > m_mcParticles
Store array.
TTree * m_treeKLM
tree for klm data
Float_t m_ECLDeltaTime
KlId for that object.
Float_t m_ECLPhi
measured phi
Float_t m_KLMhitDepth
hit depth in KLM, distance to IP
Float_t m_ECLZ51
Zernike moment 5,1 see Belle2 note on that.
Float_t m_KLMTrackClusterSepAngle
angle between trach momentum and cluster (measured from ip)
Float_t m_KLMPhi
measured phi
virtual ~DataWriterModule()
Destructor.
Float_t m_ECLMCMom
MC particle momentum; -999 if not MCparticle.
Float_t m_KLMavInterClusterDist
average distance between all KLM clusters
Float_t m_KLMECLEerror
uncertainty on E in associated ECL cluster
Float_t m_KLMMCMom
momentum of matched mc particle
Float_t m_KLMKLid
KlId for that object.
Float_t m_KLMenergy
Energy deposit in KLM (0.2 GeV * nHitCells)
Float_t m_KLMMCPrimaryPDG
pdg code of MCparticles mother, for example pi0 for some gammas
Float_t m_KLMECLZ40
zernike moment 4,0 of closest ecl cluster
Float_t m_KLMeclFlag
ecl flag for belle comparision
Float_t m_KLMECLUncertaintyPhi
phi uncertainty oof closeest ecl cluster
Float_t m_KLMECLE
energy measured in associated ECL cluster
Float_t m_ECLMCTheta
MC particle momentum; -999 if not MCparticle.
Float_t m_KLMtrackFlag
track flag for belle comparision
ECL cluster data.
Definition ECLCluster.h:27
double getDeltaL() const
Return deltaL.
Definition ECLCluster.h:265
double getUncertaintyTheta() const
Return Uncertainty on Theta of Shower.
Definition ECLCluster.h:325
double getE9oE21() const
Return E9/E21 (shower shape variable).
Definition ECLCluster.h:280
double getEnergy(EHypothesisBit hypothesis) const
Return Energy (GeV).
Definition ECLCluster.cc:23
double getUncertaintyEnergy() const
Return Uncertainty on Energy of Shower.
Definition ECLCluster.h:322
double getMinTrkDistance() const
Get distance between cluster COG and track extrapolation to ECL.
Definition ECLCluster.h:259
double getZernikeMVA() const
Return MVA based hadron/photon value based on Zernike moments (shower shape variable).
Definition ECLCluster.h:274
double getUncertaintyPhi() const
Return Uncertainty on Phi of Shower.
Definition ECLCluster.h:328
double getAbsZernike40() const
Return Zernike moment 40 (shower shape variable).
Definition ECLCluster.h:268
unsigned short getHypotheses() const
Return hypothesis (expert only, this returns a bti pattern).
Definition ECLCluster.h:244
double getDeltaTime99() const
Return cluster delta time 99.
Definition ECLCluster.h:301
double getTime() const
Return cluster time.
Definition ECLCluster.h:298
EHypothesisBit
The hypothesis bits for this ECLCluster (Connected region (CR) is split using this hypothesis.
Definition ECLCluster.h:31
@ c_nPhotons
CR is split into n photons (N1)
Definition ECLCluster.h:41
@ c_neutralHadron
CR is reconstructed as a neutral hadron (N2)
Definition ECLCluster.h:43
double getAbsZernike51() const
Return Zernike moment 51 (shower shape variable).
Definition ECLCluster.h:271
KLM cluster data.
Definition KLMCluster.h:29
Klong identification (KlId) datastore object to store results from KlId calculations.
Definition KlId.h:23
double getKlId() const
get the klong classifier output
Definition KlId.cc:29
A Class to store the Monte Carlo particle information.
Definition MCParticle.h:32
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
Module()
Constructor.
Definition Module.cc:30
Store one Track-KLMCluster separation as a ROOT object.
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Helper functions for all klid modules to improve readability of the code.
Definition KlId.h:27
bool isKLMClusterSignal(const Belle2::KLMCluster &cluster, float mcWeigthCut=0.66)
checks if a cluster is signal under the mcWeightcondition (mcWeight = energy deposition)
Definition KlId.h:100
int mcParticleIsBeamBKG(const Belle2::MCParticle *part)
return if MCparticle is beambkg
Definition KlId.h:68
std::tuple< const Belle2::KLMCluster *, double, double > findClosestKLMCluster(const ROOT::Math::XYZVector &klmClusterPosition)
find nearest KLMCluster, tis distance and the av intercluster distance
Definition KlId.h:236
bool isECLClusterSignal(const Belle2::ECLCluster &cluster, float mcWeigthCut=0.66)
checks if a cluster is signal under the mcWeightcondition (mcWeight = energy deposition)
Definition KlId.h:114
std::pair< Belle2::ECLCluster *, double > findClosestECLCluster(const ROOT::Math::XYZVector &klmClusterPosition, const Belle2::ECLCluster::EHypothesisBit eclhypothesis=Belle2::ECLCluster::EHypothesisBit::c_neutralHadron)
Find the closest ECLCluster with a neutral hadron hypothesis, and return it with its distance.
Definition KlId.h:203
int getPrimaryPDG(Belle2::MCParticle *part)
return if mc particles primary pdg.
Definition KlId.h:151
int mcParticleIsKlong(Belle2::MCParticle *part)
return the mc hirachy of the klong 0:not a klong 1:final particle, 2: klong is mother etc
Definition KlId.h:78
Abstract base class for different kinds of events.
STL namespace.