14 #include <framework/datastore/StoreArray.h>
15 #include <mdst/dataobjects/MCParticle.h>
16 #include <mdst/dataobjects/KLMCluster.h>
17 #include <mdst/dataobjects/ECLCluster.h>
18 #include <mdst/dataobjects/TrackFitResult.h>
21 #include <framework/logging/Logger.h>
22 #include <tracking/dataobjects/RecoTrack.h>
23 #include <genfit/Exception.h>
33 const TVector3& pos = cluster.getClusterPosition();
37 const TVector3& trackPos = track.getPosition();
39 if (trackPos.Angle(pos) < angle) {
40 B2DEBUG(20,
"BelleFlagTracklAngle::" << trackPos.Angle(pos));
51 const TVector3& pos = cluster.getClusterPosition();
56 const TVector3& clusterPos = eclcluster.getClusterPosition();
58 if (clusterPos.Angle(pos) < angle) {
59 B2DEBUG(20,
"BelleFlagECLAngle::" << clusterPos.Angle(pos));
69 if (part ==
nullptr) {
79 unsigned int hirachy_counter = 0;
86 if (part -> getPDG() == 130) {
87 return hirachy_counter;
89 if ((part -> getMother() ==
nullptr)) {
92 part = part -> getMother();
101 const auto mcParticleWeightPair = cluster.getRelatedToWithWeight<
Belle2::MCParticle>();
103 if (!part) {
return false; }
104 float mcWeight = mcParticleWeightPair.second;
115 const auto mcParticleWeightPair = cluster.getRelatedToWithWeight<
Belle2::MCParticle>();
117 if (!part) {
return false; }
118 float mcWeight = mcParticleWeightPair.second;
131 if (std::abs(part -> getPDG()) == pdg) {
134 if ((part -> getMother() ==
nullptr)) {
137 part = part -> getMother();
189 if ((part -> getMother() ==
nullptr)) {
192 part = part -> getMother();
195 return part ->getPDG();
207 double closestECLAngleDist = 1e10;
212 int indexOfClosestCluster = -1;
215 if (eclcluster.hasHypothesis(eclhypothesis)) {
217 const TVector3& eclclusterPos = eclcluster.getClusterPosition();
218 double angularDist = eclclusterPos.Angle(klmClusterPosition);
219 if (angularDist < closestECLAngleDist) {
220 closestECLAngleDist = angularDist;
221 indexOfClosestCluster = index;
226 if (indexOfClosestCluster > -1)
227 closestECL = eclClusters[indexOfClosestCluster];
229 return std::make_pair(closestECL, closestECLAngleDist);
234 std::tuple<const Belle2::KLMCluster*, double, double>
findClosestKLMCluster(
const TVector3& klmClusterPosition)
239 double closestKLMDist = 1e10;
240 double avInterClusterDist = 0;
241 double nKLMCluster = klmClusters.
getEntries();
243 if (nKLMCluster > 1) {
245 unsigned int index = 0;
246 unsigned int indexOfClosestCluster = 0;
249 const TVector3& nextClusterPos = nextCluster.getClusterPosition();
250 const TVector3& clustDistanceVec = nextClusterPos - klmClusterPosition;
252 double nextClusterDist = clustDistanceVec.Mag2();
253 avInterClusterDist = avInterClusterDist + nextClusterDist;
255 if ((nextClusterDist < closestKLMDist) and not(nextClusterDist == 0)) {
256 closestKLMDist = nextClusterDist ;
257 indexOfClosestCluster = index;
262 closestKLM = klmClusters[indexOfClosestCluster];
263 avInterClusterDist = avInterClusterDist / (1. * nKLMCluster);
266 return std::make_tuple(closestKLM, closestKLMDist, avInterClusterDist);
271 std::tuple<Belle2::RecoTrack*, double, std::unique_ptr<const TVector3> >
findClosestTrack(
const TVector3& clusterPosition,
275 double oldDistance = 10000000;
277 TVector3
poca = TVector3(0, 0, 0);
285 state_for_cut = track.getMeasuredStateOnPlaneFromFirstHit();
288 if (clusterPosition.Angle(state_for_cut.getPos()) < cutAngle) {
289 state = track.getMeasuredStateOnPlaneFromLastHit();
290 state.extrapolateToPoint(clusterPosition);
291 const TVector3& trackPos = state.getPos();
293 const TVector3& distanceVecCluster = clusterPosition - trackPos;
294 double newDistance = distanceVecCluster.Mag2();
297 if (newDistance < oldDistance) {
298 oldDistance = newDistance;
299 closestTrack = &track;
307 if (not closestTrack) {
308 return std::make_tuple(closestTrack, oldDistance, std::unique_ptr<const TVector3>(
nullptr));
310 return std::make_tuple(closestTrack, oldDistance, std::unique_ptr<const TVector3>(
new TVector3(poca)));