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