Belle II Software  release-06-00-14
GFTC2SPTCConverterModule.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 <tracking/modules/spacePointCreator/GFTC2SPTCConverterModule.h>
10 #include <tracking/spacePointCreation/SpacePointTrackCand.h>
11 
12 #include <framework/dataobjects/EventMetaData.h>
13 #include <framework/datastore/StoreObjPtr.h>
14 #include <framework/gearbox/Const.h>
15 
16 #include <algorithm> // count, sort, etc...
17 
18 #include <boost/tuple/tuple_comparison.hpp>
19 
20 #include <pxd/dataobjects/PXDTrueHit.h>
21 #include <svd/dataobjects/SVDTrueHit.h>
22 
23 using namespace std;
24 using namespace Belle2;
25 
26 REG_MODULE(GFTC2SPTCConverter)
27 
29  Module()
30 {
31  setDescription("Module for converting genfit::TrackCands (e.g. from TrackFinderMCTruth) to SpacePointTrackCands.");
32  setPropertyFlags(c_ParallelProcessingCertified);
33 
34  addParam("PXDClusters", m_PXDClusterName, "PXDCluster collection name", string(""));
35  addParam("SVDClusters", m_SVDClusterName, "SVDCluster collection name", string(""));
36  addParam("genfitTCName", m_genfitTCName, "Name of container of genfit::TrackCands", string(""));
37  addParam("SpacePointTCName", m_SPTCName,
38  "Name of the container under which SpacePointTrackCands will be stored in the DataStore (NOTE: These SpaceTrackCands are not checked for curling behaviour, but are simply converted and stored!)",
39  string(""));
40 
41  addParam("SingleClusterSVDSP", m_SingleClusterSVDSPName,
42  "Single Cluster SVD SpacePoints collection name. NOTE: This StoreArray will be searched for SpacePoints only if 'useSingleClusterSP' is set to true!",
43  string("SVDSpacePoints"));
44  addParam("NoSingleClusterSVDSP", m_NoSingleClusterSVDSPName,
45  "Non Single Cluster SVD SpacePoints collection name. This StoreArray will be searched for SpacePoints", string("SVDSpacePoints"));
46  addParam("PXDClusterSP", m_PXDClusterSPName, "PXD Cluster SpacePoints collection name.", string("PXDSpacePoints"));
47 
48  addParam("minNDF", m_PARAMminNDF,
49  "Minimum number of degrees of freedom a SpacePointTrackCand has to contain in order to get registered in the DataStore. If set to 0, any number is accepted",
50  0);
51 
52  addParam("checkTrueHits", m_PARAMcheckTrueHits,
53  "Set to true if you want TrueHits of Clusters forming a SpacePoint (e.g. SVD) to be checked for equality", false);
54  addParam("useSingleClusterSP", m_PARAMuseSingleClusterSP,
55  "Set to true if you want to use singleCluster SVD SpacePoints if no doubleCluster SVD SpacePoint can be found. NOTE: this gets overriden if 'skipCluster' is set to true!",
56  true);
57  addParam("checkNoSingleSVDSP", m_PARAMcheckNoSingleSVDSP,
58  "Set to false if you want to disable the initial check for the StoreArray of Non Single Cluster SVD SpacePoints. NOTE: The module will still search for these SpacePoints first, so you have to make sure you are not registering SpacePoints under the StoreArray with the NoSingleClusterSVDSP name! (Disable the module that registers these SpacePoints)",
59  true);
60  addParam("skipCluster", m_PARAMskipCluster,
61  "Set to true if you only want to skip the Clusters for which no appropriate SpacePoints can be found, instead of aborting the conversion of the whole GFTC when such a case occurs. NOTE: setting this to true automatically sets 'useSingleClusterSP' to false!",
62  false);
63 
64  initializeCounters(); // NOTE: they get initialized in initialize again!!
65 }
66 
67 // ------------------------------ INITIALIZE ---------------------------------------
68 void GFTC2SPTCConverterModule::initialize()
69 {
70  B2INFO("GFTC2SPTCConverter -------------- initialize() ---------------------");
71  // initialize Counters
72  initializeCounters();
73 
74  // check if all required StoreArrays are here
75  StoreArray<PXDCluster> PXDClusters(m_PXDClusterName); PXDClusters.isRequired(m_PXDClusterName);
76  StoreArray<SVDCluster> SVDClusters(m_SVDClusterName); SVDClusters.isRequired(m_SVDClusterName);
77  if (m_PARAMuseSingleClusterSP) {
78  StoreArray<SpacePoint> SCSPs(m_SingleClusterSVDSPName);
79  SCSPs.isRequired(m_SingleClusterSVDSPName);
80  }
81  if (m_PARAMcheckNoSingleSVDSP) {
82  StoreArray<SpacePoint> nSCSPs(m_NoSingleClusterSVDSPName);
83  nSCSPs.isRequired(m_NoSingleClusterSVDSPName);
84  }
85  StoreArray<SpacePoint> pxdSPs(m_PXDClusterSPName);
86  pxdSPs.isRequired(m_PXDClusterSPName);
87 
88  StoreArray<genfit::TrackCand> gfTrackCand(m_genfitTCName);
89  gfTrackCand.isRequired(m_genfitTCName);
90 
91  // registering StoreArray for SpacePointTrackCand
92  StoreArray<SpacePointTrackCand> spTrackCand(m_SPTCName);
93  spTrackCand.registerInDataStore(m_SPTCName, DataStore::c_ErrorIfAlreadyRegistered);
94 
95  // register Relation to genfit::TrackCand
96  spTrackCand.registerRelationTo(gfTrackCand);
97 
98  // CAUTION: if the StoreArray of the TrueHits is named, this check fails!!!
99  if (m_PARAMcheckTrueHits) {
100  StoreArray<PXDTrueHit> PXDTrueHits; PXDTrueHits.isRequired();
101  StoreArray<SVDTrueHit> SVDTrueHits; SVDTrueHits.isRequired();
102  }
103 
104  if (m_PARAMminNDF < 0) {
105  B2WARNING("'minNDF' is set to a value below 0. Resetting to 0!");
106  m_PARAMminNDF = 0;
107  }
108 
109 }
110 
111 // ------------------------------------- EVENT -------------------------------------------------------
112 void GFTC2SPTCConverterModule::event()
113 {
114  StoreObjPtr<EventMetaData> eventMetaDataPtr("EventMetaData", DataStore::c_Event);
115  const int eventCounter = eventMetaDataPtr->getEvent();
116  B2DEBUG(10, "GFTC2SPTCConverter::event(). Processing event " << eventCounter << " --------");
117 
118  StoreArray<genfit::TrackCand> mcTrackCands(m_genfitTCName);
119  StoreArray<SpacePointTrackCand> spacePointTrackCands(m_SPTCName); // output StoreArray
120 
121  StoreArray<PXDCluster> pxdClusters(m_PXDClusterName);
122  StoreArray<SVDCluster> svdClusters(m_SVDClusterName);
123 
124  int nTCs = mcTrackCands.getEntries();
125 
126  B2DEBUG(15, "Found " << nTCs << " genfit::TrackCands in StoreArray " << mcTrackCands.getName());
127 
128  for (int iTC = 0; iTC < nTCs; ++iTC) {
129  genfit::TrackCand* trackCand = mcTrackCands[iTC];
130  m_genfitTCCtr += 1;
131 
132  B2DEBUG(10,
133  "================================================================================\nNow processing genfit::TrackCand "
134  << iTC << ".");
135  if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 15, PACKAGENAME())) { trackCand->Print(); } // prints to stdout
136  try {
137  // get the converted SPTC
138  std::pair<const SpacePointTrackCand, conversionStatus> spacePointTC = createSpacePointTC(trackCand, pxdClusters, svdClusters);
139  // check if the SPTC contains enough SpacePoints
140  if (spacePointTC.second == 0) {
141  if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 50, PACKAGENAME())) spacePointTC.first.print(50);
142  SpacePointTrackCand* newSPTC = spacePointTrackCands.appendNew(spacePointTC.first);
143  m_SpacePointTCCtr++;
144  newSPTC->addRelationTo(trackCand);
145  B2DEBUG(10, "Added new SpacePointTrackCand to StoreArray " << spacePointTrackCands.getName());
146  } else {
147  B2DEBUG(10, "The conversion failed due to: " << spacePointTC.second);
148  increaseFailCounter(spacePointTC.second);
149  }
150  } catch (std::runtime_error& anE) { // catch all Belle2 exceptions with this!
151  B2ERROR("Caught exception in creation of SpacePointTrackCand: " << anE.what());
152  } catch (...) { // catch the rest
153  B2ERROR("Caught undefined exception in creation of SpacePointTrackCand!"); // COULDDO: rethrow exception as something is fishy if this happens
154  }
155  }
156 }
157 
158 // -------------------------------- TERMINATE --------------------------------------------------------
159 void GFTC2SPTCConverterModule::terminate()
160 {
161  stringstream generalOutput;
162  generalOutput << "GFTC2SPTCConverter::terminate(): got " << m_genfitTCCtr << " GFTCs and created " << m_SpacePointTCCtr <<
163  " SPTCs. ";
164  if (m_abortedLowNDFCtr) generalOutput << "For " << m_abortedLowNDFCtr << " SPTCs the NDF was below " << m_PARAMminNDF << "\n";
165  // NEW output
166  B2INFO(generalOutput.str());
167  if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 1, PACKAGENAME())) {
168  stringstream verboseOutput;
169  verboseOutput << "counter variables: ";
170  verboseOutput << "abortedNoSP: " << m_abortedNoSPCtr << ", abortedUnsuitableGFTC: " << m_abortedUnsuitableTCCtr <<
171  ", abortedNoValidSP: " << m_abortedNoValidSPCtr;
172  if (m_PARAMcheckTrueHits) verboseOutput << ", abortedTrueHit: " << m_abortedTrueHitCtr;
173  if (m_PARAMskipCluster) {
174  verboseOutput << ", skippedCluster: " << m_skippedCluster << ", skippedPXDnoSP: " << m_skippedPXDnoSPCtr;
175  verboseOutput << ", skippedSVDnoSP: " << m_skippedSVDnoSPCtr << ", skippedPXDnoTH: " << m_skippedPXDnoTHCtr << ", skippedSVDnoTH: "
176  << m_skippedSVDnoTHCtr;
177  verboseOutput << ", skippedPXDunsuitable: " << m_skippedPXDunsuitableCtr << ", skippedSVDunsuitable: " << m_skippedSVDunsuitableCtr;
178  verboseOutput << ", skippedPXDnoValidSP " << m_skippedPXDnoValidSPCtr << ", skippedSVDnoValidSP " << m_skippedSVDnoValidSPCtr;
179  }
180  if (m_nonSingleSPCtr) verboseOutput << ", nonSingleSP " << m_nonSingleSPCtr;
181  verboseOutput << ", noTwoClusterSP: " << m_noTwoClusterSPCtr << ", singleClusterSVDSP: " << m_singleClusterSPCtr;
182  B2DEBUG(1, verboseOutput.str());
183  }
184  if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 2, PACKAGENAME())) {
185  stringstream explanation;
186  explanation << "explanation of counter variables (key words only):\n";
187  explanation << "NoSP -> Found no related SpacePoint to a Cluster\n";
188  explanation << "Unsuitable -> Cluster combination of SpacePoint was not in consecutive order in GFTC\n";
189  explanation << "NoValidSP -> Cluster combination of SpacePoint was not contained in GFTC\n";
190  if (m_PARAMcheckTrueHits) explanation << "TrueHit/noTH -> found no related TrueHit to a SpacePoint\n";
191  if (m_nonSingleSPCtr) explanation << "nonSingleSP -> more than one singleCluster SpacePoint related to a Cluster\n";
192  explanation << "noTwoClusterSP -> found no two Cluster SpacePoint\n";
193  explanation << "singleClusterSVDSP -> number of tries to add a singleCluster SpacePoint for latter cases\n";
194  B2DEBUG(2, explanation.str());
195  }
196 }
197 
198 // ================================================================= CREATE SPACEOINT TRACKCAND =========================================================================================
199 std::pair<const Belle2::SpacePointTrackCand, GFTC2SPTCConverterModule::conversionStatus>
200 GFTC2SPTCConverterModule::createSpacePointTC(const genfit::TrackCand* genfitTC, const StoreArray<PXDCluster>& pxdClusters,
201  const StoreArray<SVDCluster>& svdClusters)
202 {
203  m_NDF = 0; // reset for every trackCand
204  std::vector<HitInfo<SpacePoint> > tcSpacePoints;
205  // cppcheck-suppress unreadVariable
206  conversionStatus convStatus = c_noFail; // part of return
207 
208  int nHits = genfitTC->getNHits();
209  B2DEBUG(15, "genfit::TrackCand contains " << nHits << " hits");
210 
211  // for easier handling fill a taggedPair (typedef) to distinguish whether a hit has already been used or if is still open for procession
212  std::vector<flaggedPair<int> > fHitIDs;
213  for (int i = 0; i < nHits; ++i) {
214  auto aHit = genfitTC->getHit(i);
215  flaggedPair<int> aPair(false, aHit->getDetId(), aHit->getHitId());
216  fHitIDs.push_back(aPair);
217  }
218 
219  bool usedSingleCluster = false;
220  // now loop over all hits and add them appropriately
221  for (int iTCHit = 0; iTCHit < nHits; ++iTCHit) {
222  genfit::TrackCandHit* aTCHit = genfitTC->getHit(iTCHit);
223  double sortingParam = aTCHit->getSortingParameter();
224 
225  B2DEBUG(20, "Processing TrackCandHit " << iTCHit << " of " << nHits);
226  if (fHitIDs[iTCHit].get<0>()) { // check if this hit has already been used, if not process
227  B2DEBUG(60, "This hit has already been added to the SpacePointTrackCand via a SpacePoint and will not be processed again");
228  } else {
229  std::pair<SpacePoint*, conversionStatus> aSpacePoint = processTrackCandHit(aTCHit, pxdClusters, svdClusters, fHitIDs, iTCHit);
230 
231  // if there is more than one single Cluster SpacePoint related to a Cluster, add one
232  if (aSpacePoint.first != nullptr && (aSpacePoint.second >= 0 || aSpacePoint.second == c_nonSingleSP)) {
233  if (aSpacePoint.second == c_singleClusterSP) usedSingleCluster = true;
234  tcSpacePoints.push_back({sortingParam, aSpacePoint.first});
235  B2DEBUG(60, "Added SpacePoint " << aSpacePoint.first->getArrayIndex() << " from Array " << aSpacePoint.first->getArrayName() <<
236  " to tcSpacePoints");
237  m_NDF += getNDF(aSpacePoint.first);
238  } else {
239  convStatus = aSpacePoint.second;
240  B2DEBUG(60, "There was an error during conversion: for Hit " << iTCHit << ": " << convStatus);
241  if (!m_PARAMskipCluster) {
242  B2DEBUG(10,
243  "There was an error during conversion for a GFTC. 'skipCluster' is set to false, hence this trackCand will not be converted!");
244 
245  // return empty SPTC if there is a conversion error and splitCurlers is set to false
246  return std::make_pair(SpacePointTrackCand(), convStatus);
247  }
248  }
249  }
250  }
251 
252  B2DEBUG(20, "NDF for this SpacePointTrackCand: " << m_NDF);
253  if (m_NDF < m_PARAMminNDF) {
254  B2DEBUG(10, "The created SpacePointTrackCand has not enough NDF: NDF is " << m_NDF << " but 'minNDF' is set to " << m_PARAMminNDF);
255  return std::make_pair(SpacePointTrackCand(), c_lowNDF);
256  }
257 
258  // check if all hits have been used
259  bool usedAllHits = checkUsedAllHits(fHitIDs);
260  if (!usedAllHits && !m_PARAMskipCluster) {
261  B2WARNING("There is at least one TrackCandHit that has not been marked as used although 'skipCluster' is set to false"); // write warning here, because if this happens something has gone wrong
262  throw UnusedHits();
263  }
264 
265  // create a vector of SpacePoint* and one with sorting Parameters to add to the SpacePointTrackCand
266  std::vector<const SpacePoint*> spacePoints;
267  std::vector<double> sortingParams;
268  for (const HitInfo<SpacePoint>& aSP : tcSpacePoints) {
269  spacePoints.push_back(aSP.second);
270  sortingParams.push_back(aSP.first);
271  }
272 
273  SpacePointTrackCand spacePointTC = SpacePointTrackCand(spacePoints, genfitTC->getPdgCode(), genfitTC->getChargeSeed(),
274  genfitTC->getMcTrackId());
275  spacePointTC.set6DSeed(genfitTC->getStateSeed());
276  spacePointTC.setCovSeed(genfitTC->getCovSeed());
277  spacePointTC.setSortingParameters(sortingParams);
278 
279  if (m_PARAMcheckTrueHits) { spacePointTC.addRefereeStatus(SpacePointTrackCand::c_checkedTrueHits); }
280  if (!usedAllHits) { spacePointTC.addRefereeStatus(SpacePointTrackCand::c_omittedClusters); }
281  if (usedSingleCluster) { spacePointTC.addRefereeStatus(SpacePointTrackCand::c_singleClustersSPs); }
282 
283  return std::make_pair(spacePointTC, c_noFail);
284 }
285 
286 // ============================================================================== PROCESS TRACKCANDHIT ===============================================================================================
287 std::pair<Belle2::SpacePoint*, GFTC2SPTCConverterModule::conversionStatus>
288 GFTC2SPTCConverterModule::processTrackCandHit(genfit::TrackCandHit* hit, const StoreArray<PXDCluster>& pxdClusters,
289  const StoreArray<SVDCluster>& svdClusters,
290  std::vector<flaggedPair<int> >& flaggedHitIDs, int iHit)
291 {
292  int detID = hit->getDetId();
293  int hitID = hit->getHitId();
294  int planeID = hit->getPlaneId(); // not used at the moment (except for debug output)
295  B2DEBUG(60, "Processing TrackCandHit " << iHit << " with detID: " << detID << ", hitID: " << hitID << ", planeID: " << planeID);
296 
297  std::pair<SpacePoint*, conversionStatus> returnSP = { nullptr, c_noFail }; // default return, optimistically assuming no fail
298 
299  if (detID == Const::PXD) {
300  const PXDCluster* aCluster = pxdClusters[hitID];
301  returnSP = getSpacePoint<PXDCluster, PXDTrueHit>(aCluster, flaggedHitIDs, iHit, true, m_PXDClusterSPName);
302  if (m_PARAMskipCluster) { increaseSkippedCounter(returnSP.second, aCluster); } // have to do this here at the moment -> COULDDO; change increaseSkippedCounter, such that it can be used without a Cluster
303  if (returnSP.second == c_singleClusterSP) returnSP.second = c_noFail; // getSpacePoint returns singleClusterSP for PXDs!
304  } else if (detID == Const::SVD) {
305  const SVDCluster* aCluster = svdClusters[hitID];
306  returnSP = getSpacePoint<SVDCluster, SVDTrueHit>(aCluster, flaggedHitIDs, iHit, false, m_NoSingleClusterSVDSPName);
307  if (m_PARAMskipCluster) { increaseSkippedCounter(returnSP.second, aCluster); }
308  } else {
309  throw SpacePointTrackCand::UnsupportedDetType();
310  }
311  return returnSP;
312 }
313 
314 // ============================================================================ GET SPACEPOINT =======================================================================================================
315 template<typename ClusterType, typename TrueHitType>
316 std::pair<Belle2::SpacePoint*, GFTC2SPTCConverterModule::conversionStatus>
317 GFTC2SPTCConverterModule::getSpacePoint(const ClusterType* cluster, std::vector<flaggedPair<int> >& flaggedHitIDs, int iHit,
318  bool singleCluster, std::string arrayName)
319 {
320  std::pair<SpacePoint*, conversionStatus> spacePoint = {nullptr, c_noFail}; // default return. be optimistic and assume that there are no problems!
321 
322  B2DEBUG(70, "Trying to find a related SpacePoint in StoreArray " << arrayName << " for Cluster " << cluster->getArrayIndex() <<
323  " from Array " << cluster->getArrayName());
324  RelationVector<SpacePoint> spacePoints = cluster->template getRelationsFrom<SpacePoint>(arrayName);
325  B2DEBUG(80, "Found " << spacePoints.size() << " related SpacePoints for Cluster " << cluster->getArrayIndex() << " from Array " <<
326  cluster->getArrayName()); // NOTE: .size == 0 handled later!
327 
328  if (singleCluster) { // if it is a singleCluster SpacePoint process different, then if it is a double Cluster SpacePoint
329  if (spacePoints.size() == 0) {
330  B2DEBUG(80, "Found no related (single Cluster) SpacePoint!");
331  spacePoint.second = c_foundNoSpacePoint;
332  return spacePoint;
333  }
334  if (spacePoints.size() > 1) {
335  B2ERROR("More than one single Cluster SpacePoint related to a Cluster! Returning only the first in RelationVector!");
336  spacePoint.second = c_nonSingleSP; // WARNING: returning first SpacePoint in RelationVector this way
337  m_nonSingleSPCtr++;
338  }
339  spacePoint.first = spacePoints[0];
340  // cppcheck-suppress redundantAssignment
341  spacePoint.second = c_singleClusterSP;
342  markHitAsUsed(flaggedHitIDs, iHit); // mark hit as used
343  } else {
344  // try to get the appropriate double Cluster SVD SpacePoint
345  spacePoint = findAppropriateSpacePoint<ClusterType>(spacePoints, flaggedHitIDs);
346  if (spacePoint.first == nullptr) {
347  B2DEBUG(80, "Did not find an appropriate double Cluster SpacePoint for Cluster " << cluster->getArrayIndex() << " from Array " <<
348  cluster->getArrayName() << ". Reason for failure: " << spacePoint.second);
349  m_noTwoClusterSPCtr++;
350  if (m_PARAMuseSingleClusterSP) {
351  B2DEBUG(80, "Trying to get a single Cluster SpacePoint now!");
352  m_singleClusterSPCtr++;
353  // get single Cluster SVD PacePoint if desired. WARNING: hardcoded to SVD here at the moment!
354  return getSpacePoint<ClusterType, TrueHitType>(cluster, flaggedHitIDs, iHit, true, m_SingleClusterSVDSPName);
355  }
356  }
357  }
358 
359  if (m_PARAMcheckTrueHits && spacePoint.first != nullptr) { // only do the TrueHit check if there is actually something to check
360  if (!foundRelatedTrueHit<TrueHitType>(spacePoint.first)) { spacePoint.second = c_foundNoTrueHit; }
361  }
362 
363  return spacePoint;
364 }
365 
366 // ============================================================= FIND APPROPRIATE SPACEPOINT ===============================================================================
367 template<typename ClusterType>
368 std::pair<Belle2::SpacePoint*, GFTC2SPTCConverterModule::conversionStatus>
369 GFTC2SPTCConverterModule::findAppropriateSpacePoint(const Belle2::RelationVector<Belle2::SpacePoint>& spacePoints,
370  std::vector<flaggedPair<int> >& flaggedHitIDs)
371 {
372  std::pair<SpacePoint*, conversionStatus> returnSP = { nullptr, c_noFail }; // default return, be optimistc and assume that there are no problems
373  B2DEBUG(100, "Trying to find an appropriate SpacePoint from RelationVector with " << spacePoints.size() << " entries!");
374  if (spacePoints.size() == 0) {
375  B2DEBUG(80, "There are no spacePoints to choose of!");
376  returnSP.second = c_foundNoSpacePoint;
377  return returnSP;
378  }
379 
380  // WARNING: TEL cluster will do something wrong here!
381  int detID = spacePoints[0]->getType() == VXD::SensorInfoBase::SVD ? Const::SVD : Const::PXD;
382 
383  // .first: Cluster Combination exists in the TrackCand, .second: Cluster Combination exists and has not yet been used
384  std::vector<std::pair<bool, bool> > existAndValidSP;
385  // .first: index of SpacePoint in RelationVector, .second: Index of Cluster in the genfit::TrackCand
386  std::vector<std::pair<int, int> > clusterPositions;
387 
388  // loop over all SpacePoints to look if their Cluster combination is vali (or existing but used, etc...)
389  for (unsigned int iSP = 0; iSP < spacePoints.size(); ++iSP) {
390  const SpacePoint* aSP = spacePoints[iSP];
391  B2DEBUG(100, "Processing SpacePoint " << iSP + 1 << " of " << spacePoints.size() << " with Index " << aSP->getArrayIndex() <<
392  " in StoreArray " << aSP->getArrayName());
393 
394  bool bothValid = true; // assume true, change to false if the second Cluster cannot be found in the genfit::TrackCand
395  bool foundBoth = true; // assume true, change to false if the second Cluster can be found, but is already used
396 
397  std::vector<int> clusterInds = getClusterIndices<SVDCluster>(aSP, m_SVDClusterName);
398  for (int index : clusterInds) {
399  // get the valid and existing position for cluster. TODO: rename variable!
400  std::pair<int, int> existAndValidClPos = checkExistAndValid(index, detID, flaggedHitIDs);
401  if (existAndValidClPos.first < 0) { // check if cluster is valid and/or exists
402  bothValid = false;
403  if (existAndValidClPos.second < 0) foundBoth = false;
404  }
405 
406  // push_back the index in the relationVector and the found valid position
407  clusterPositions.push_back({iSP, existAndValidClPos.first});
408  B2DEBUG(999, "clusterInd: " << index << " checkExistAndValid.first: " << existAndValidClPos.first << ", .second: " <<
409  existAndValidClPos.second << " bothValid/foundBoth: " << bothValid << "/" << foundBoth);
410  }
411  // push_back the determined values
412  existAndValidSP.push_back({foundBoth, bothValid});
413  B2DEBUG(100, "Cluster combination of SpacePoint " << aSP->getArrayIndex() << " is contained in genfit::TrackCand: " << foundBoth <<
414  ". SpacePoint is valid: " << bothValid);
415  }
416 
417  int relVecPosition = getAppropriateSpacePointIndex(existAndValidSP, clusterPositions);
418  if (relVecPosition < 0) {
419  returnSP.second = getFailEnum(relVecPosition);
420  return returnSP;
421  }
422 
423  B2DEBUG(100, "SpacePoint " << spacePoints[relVecPosition]->getArrayIndex() <<
424  " is the appropriate SpacePoint of all checked SpacePoints! The positions inside the GFTC are: " << clusterPositions.at(
425  relVecPosition * 2).second << " and " << clusterPositions.at(relVecPosition * 2 + 1).second);
426  markHitAsUsed(flaggedHitIDs, clusterPositions.at(relVecPosition * 2).second);
427  markHitAsUsed(flaggedHitIDs, clusterPositions.at(relVecPosition * 2 + 1).second);
428  returnSP.first = spacePoints[relVecPosition];
429 
430  return returnSP;
431 }
432 
433 // ========================================================================= GET CLUSTER INDICES ===========================================================================
434 template<typename ClusterType>
435 std::vector<int> GFTC2SPTCConverterModule::getClusterIndices(const Belle2::SpacePoint* spacePoint, std::string storeArrayName)
436 {
437  std::vector<int> clusterInds;
438  RelationVector<ClusterType> relClusters = spacePoint->getRelationsTo<ClusterType>(storeArrayName);
439  B2ASSERT("Too many clusters. There are " << relClusters.size() << " clusters.", relClusters.size() < 3);
440 
441  stringstream clusterStream;
442  for (const ClusterType& aCluster : relClusters) {
443  clusterInds.push_back(aCluster.getArrayIndex());
444  clusterStream << aCluster.getArrayIndex() << " ";
445  }
446 
447  B2DEBUG(499, "getClusterIndices(SpacePoint " << spacePoint->getArrayIndex() << "," << spacePoint->getArrayName() <<
448  "): clusters are: " << clusterStream.str());
449 
450  return clusterInds;
451 }
452 
453 // ========================================================================== CHECK EXIST AND VALID =============================================================================
454 std::pair<int, int> GFTC2SPTCConverterModule::checkExistAndValid(int clusterInd, int detID,
455  std::vector<flaggedPair<int> >& flaggedHitIDs)
456 {
457  B2DEBUG(499, "Now checking if Cluster " << clusterInd << " is valid");
458  std::pair<int, int> positions = { -1, -1}; // return Vector
459 
460  // flaggedPair for finding valid IDs: hit is in genfit::TrackCand and has not yet been used by another SpacePoint
461  flaggedPair<int> validID(false, detID, clusterInd);
462  flaggedPair<int> existingID(true, detID, clusterInd); // flaggedPair for finding existing but used fHitIDs
463 
464  // find the positions of these two pairs in flaggedHitIDs
465  // position in "normal array indexing" (i.e. starting at 0, ending at vector.size() -1 )
466  unsigned int validPos = std::find(flaggedHitIDs.begin(), flaggedHitIDs.end(), validID) - flaggedHitIDs.begin();
467  unsigned int existingPos = std::find(flaggedHitIDs.begin(), flaggedHitIDs.end(), existingID) - flaggedHitIDs.begin();
468 
469  B2DEBUG(100, "validID = (" << validID.get<1>() << "," << validID.get<2>() << "), found at position " << validPos <<
470  ", existingID found at position " << existingPos << " of " << flaggedHitIDs.size());
471 
472  // check if these positions are still in the TrackCand
473  if (validPos < flaggedHitIDs.size()) { positions.first = validPos; }
474  if (existingPos < flaggedHitIDs.size()) { positions.second = existingPos; }
475 
476  B2DEBUG(999, "Return values, .first: " << positions.first << ", .second: " << positions.second);
477  return positions;
478 }
479 
480 // ============================================================================ GET NDF =============================================================================================
481 int GFTC2SPTCConverterModule::getNDF(Belle2::SpacePoint* spacePoint)
482 {
483  if (spacePoint == nullptr) {
484  B2ERROR("Got nullptr pointer to determine the NDF of!");
485  return 0;
486  }
487  std::pair<bool, bool> assignedHits = spacePoint->getIfClustersAssigned();
488  // if both are assigned -> NDF of SpacePoint is 2
489  if (assignedHits.first && assignedHits.second) return 2;
490  // if ONLY one is assigned -> NDF of SpacePoint is 1 (note the use of the bitwise XOR operator)
491  if (assignedHits.first ^ assignedHits.second) return 1;
492  return 0;
493 }
494 
495 // =========================================================================== GET APPROPRIATE SPACEPOINT INDEX =====================================================================
496 int GFTC2SPTCConverterModule::getAppropriateSpacePointIndex(const std::vector<std::pair<bool, bool> >& existAndValidSPs,
497  const std::vector<std::pair<int, int> >& clusterPositions)
498 {
499  // get the number of existing but used and the number of valid SpacePoints
500  int nExistingButUsedSP = std::count(existAndValidSPs.begin(), existAndValidSPs.end(), std::make_pair(true, false));
501  int nValidSP = std::count(existAndValidSPs.begin(), existAndValidSPs.end(), std::make_pair(true, true));
502 
503  if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 999, PACKAGENAME())) { // some very verbose output
504  stringstream output;
505  output << "content of passed vector of pairs (comma separated): ";
506  for (auto entry : existAndValidSPs) { output << entry.first << "/" << entry.second << ", "; }
507  B2DEBUG(999, output.str() << "nValidSP: " << nValidSP << " nExistingButUsedSP: " << nExistingButUsedSP);
508  }
509 
510  // if there is no valid SpacePoint, but a SpacePoint with an existing but used cluster, throw, because conversion cannot be done properly then
511  if (nValidSP < 1 && nExistingButUsedSP > 0) {
512  B2DEBUG(80,
513  "There are only Cluster Combinations where one of the Clusters is already used by another SpacePoint! This genfit::TrackCand cannot be converted properly to a SpacePointTrackCand!");
514  return c_unsuitableGFTC;
515  } else if (nValidSP < 1 && nExistingButUsedSP < 1) {
516  B2DEBUG(120, "Found no valid SpacePoint and no SpacePoint with existing but used Clusters/Hits!");
517  return c_noValidSP;
518  } // if there is at least one valid SpacePoint, check the position difference and then decide if the SP can be used!
519  else if (nValidSP > 0) {
520  // 1) index of SpacePoint in RelationVector, 2) squared position difference, 3) & 4) are positions inside genfit::TrackCand (valid positions)
521  std::vector<std::pair<int, int> > positionInfos;
522 
523  for (unsigned int iSP = 0; iSP < existAndValidSPs.size(); ++iSP) {
524  if (!existAndValidSPs.at(iSP).second) continue; // if not valid continue with next SpacePoint
525  // sign doesnot matter, comparing squared values later only!
526  int posDiff = clusterPositions.at(iSP * 2).second - clusterPositions.at(iSP * 2 + 1).second;
527 
528  B2DEBUG(200, "Difference of positions of Clusters for entry " << iSP << " is " << posDiff);
529  positionInfos.push_back(std::make_pair(iSP, posDiff * posDiff));
530  }
531 
532  // sort to find the smallest difference
533  sort(positionInfos.begin(), positionInfos.end(), [](const pair<int, int>& lTuple, const pair<int, int>& rTuple) { return lTuple.second < rTuple.second; });
534 
535  if (positionInfos.at(0).second != 1) {
536  B2DEBUG(80, "The shortest squared distance between two Clusters is " << positionInfos.at(0).second <<
537  "! This would lead to wrong ordered TrackCandHits.");
538  return c_unsuitableGFTC;
539  }
540 
541  B2DEBUG(150, "SpacePoint with index " << positionInfos.at(0).first <<
542  " is the valid SpacePoint with two Clusters in consecutive order from all valid SpacePoints.");
543  return positionInfos.at(0).first;
544  }
545 
546  return c_noValidSP; // default return is negative
547 }
548 
549 // =============================================================================== CHECK USED ALL HITS ==========================================================================
550 bool GFTC2SPTCConverterModule::checkUsedAllHits(std::vector<flaggedPair<int> >& flaggedHitIDs)
551 {
552  bool usedAll = true; // defining variable here so that all hits are checked (debug output)
553  for (unsigned int i = 0; i < flaggedHitIDs.size(); ++i) {
554  flaggedPair<int> fPair = flaggedHitIDs[i];
555  B2DEBUG(200, "Hit " << i << " with (detID,hitID): (" << fPair.get<1>() << "," << fPair.get<2>() << ") has been used: " <<
556  fPair.get<0>());
557  if (!fPair.get<0>()) {
558 // return false;
559  usedAll = false;
560  }
561  }
562 // return true;
563  return usedAll;
564 }
565 
566 // ================================================== FOUND RELATED TRUEHIT =================================================================
567 template <typename TrueHitType>
568 bool GFTC2SPTCConverterModule::foundRelatedTrueHit(const Belle2::SpacePoint* spacePoint, unsigned int allowedRelations)
569 {
570  RelationVector<TrueHitType> relTrueHits = spacePoint->getRelationsTo<TrueHitType>("ALL"); // WARNING: searching in all relations
571  if (relTrueHits.size() == 0) {
572  B2DEBUG(100, "Found no TrueHit to SpacePoint " << spacePoint->getArrayIndex() << " from Array " << spacePoint->getArrayName());
573  return false;
574  }
575  B2DEBUG(100, "Found " << relTrueHits.size() << " related TrueHits for SpacePoint " << spacePoint->getArrayIndex() << " from Array "
576  << spacePoint->getArrayName());
577  return (relTrueHits.size() <= allowedRelations);
578 }
579 
580 
581 
582 void GFTC2SPTCConverterModule::markHitAsUsed(std::vector<flaggedPair<int> >& flaggedHitIDs, int hitToMark)
583 {
584  flaggedHitIDs[hitToMark].get<0>() = true;
585  flaggedPair<int> fPair = flaggedHitIDs[hitToMark];
586  B2DEBUG(150, "Marked Hit " << hitToMark << " as used. (detID,hitID) of this hit is (" << fPair.get<1>() << "," << fPair.get<2>() <<
587  ")");
588 }
Module for converting genfit::TrackCands to SpacePointTrackCands.
conversionStatus
enum for differentiating different reasons why a conversion failed negative values mean fail!
std::pair< double, const HitType * > HitInfo
container used for storing information, that is then put into the SpacePointTrackCand
boost::tuple< bool, T, T > flaggedPair
typdef, for avoiding having a vector<bool> and a vector<pair<T,T>>
Base class for Modules.
Definition: Module.h:72
The PXD Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition: PXDCluster.h:30
Class for type safe access to objects that are referred to in relations.
size_t size() const
Get number of relations.
void addRelationTo(const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
Add a relation from this object to another object (with caching).
std::string getArrayName() const
Get name of array this object is stored in, or "" if not found.
int getArrayIndex() const
Returns this object's array index (in StoreArray), or -1 if not found.
RelationVector< TO > getRelationsTo(const std::string &name="", const std::string &namedRelation="") const
Get the relations that point from this object to another store array.
The SVD Cluster class This class stores all information about reconstructed SVD clusters.
Definition: SVDCluster.h:28
Storage for (VXD) SpacePoint-based track candidates.
void setSortingParameters(const std::vector< double > &sortParams)
set the sorting parameters
void set6DSeed(const TVectorD &state6D)
set the 6D state seed
void setCovSeed(const TMatrixDSym &cov)
set the covariance matrix seed
void addRefereeStatus(unsigned short int bitmask)
add a referee status
SpacePoint typically is build from 1 PXDCluster or 1-2 SVDClusters.
Definition: SpacePoint.h:42
std::pair< bool, bool > getIfClustersAssigned() const
Returns, if u(v)-coordinate is based on cluster information.
Definition: SpacePoint.h:162
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.
const std::string & getName() const
Return name under which the object is saved in the DataStore.
bool registerInDataStore(DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut)
Register the object/array in the DataStore.
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:246
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:216
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
Definition: StoreArray.h:140
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:95
Hit object for use in TrackCand.
Definition: TrackCandHit.h:34
Track candidate – seed values and indices.
Definition: TrackCand.h:69
int getMcTrackId() const
Get the MCT track id, for MC simulations - default value = -1.
Definition: TrackCand.h:119
int getPdgCode() const
Get the PDG code.
Definition: TrackCand.h:139
void Print(const Option_t *="") const
Write the content of all private attributes to the terminal.
Definition: TrackCand.cc:202
const TMatrixDSym & getCovSeed() const
get the covariance matrix seed (6D).
Definition: TrackCand.h:131
const TVectorD & getStateSeed() const
Returns the 6D seed state; should be in global coordinates.
Definition: TrackCand.h:134
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:650
Abstract base class for different kinds of events.