Belle II Software  release-05-01-25
CellularAutomaton.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2012 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Oliver Frost *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #pragma once
11 
12 #include <tracking/trackFindingCDC/ca/AutomatonCell.h>
13 #include <tracking/trackFindingCDC/numerics/Weight.h>
14 
15 #include <tracking/trackFindingCDC/utilities/WeightedRelation.h>
16 
17 #include <framework/logging/Logger.h>
18 
19 #include <cmath>
20 
21 namespace Belle2 {
27  namespace TrackFindingCDC {
31  template<class ACellHolder>
32  class CellularAutomaton {
33 
34  private:
36  class CycleException {};
37 
38  public:
45  ACellHolder* applyTo(const std::vector<ACellHolder*>& cellHolders,
46  const std::vector<WeightedRelation<ACellHolder>>& cellHolderRelations) const
47  {
48  B2ASSERT("Expected the relations to be sorted",
49  std::is_sorted(cellHolderRelations.begin(), cellHolderRelations.end()));
50 
51  // Set all cell states to -inf and the non permanent flags to unset.
52  prepareCellFlags(cellHolders);
53 
54  for (ACellHolder* cellHolder : cellHolders) {
55  AutomatonCell& cell = cellHolder->getAutomatonCell();
56  if (cell.hasMaskedFlag()) continue;
57  if (cell.hasCycleFlag()) continue;
58 
59  if (not cell.hasAssignedFlag()) {
60  // Mark this cell as a start point of a long path since we encountered it in
61  // at the top level of the recursion.
62  cell.setStartFlag();
63  // The flag will be unset when it appears on the _to_ side of a relation.
64  }
65 
66  try {
67  // Assignes flags and the cell state
68  getFinalCellState(cellHolder, cellHolderRelations);
69  } catch (CycleException) {
70  // TODO: Come up with some handeling for cycles.
71  // For now we continue to look for long paths in the graph
72  // hoping to find a part that does not enter the cycle.
73 
74  // Thoughts:
75  // If there is a single cycle in the graph we might be able to break it at some point.
76  // However, if there are multiple cycles intertwined things get very tricky.
77  // But can we actually distinguish the two situations?
78  }
79  }
80 
81  auto lessStartCellState = [](ACellHolder * lhs, ACellHolder * rhs) {
82  AutomatonCell& lhsCell = lhs->getAutomatonCell();
83  AutomatonCell& rhsCell = rhs->getAutomatonCell();
84  return (std::make_pair(lhsCell.hasStartFlag(), lhsCell.getCellState()) <
85  std::make_pair(rhsCell.hasStartFlag(), rhsCell.getCellState()));
86  };
87 
88  auto itStartCellHolder =
89  std::max_element(cellHolders.begin(), cellHolders.end(), lessStartCellState);
90  if (itStartCellHolder == cellHolders.end()) return nullptr;
91  if (not(*itStartCellHolder)->getAutomatonCell().hasStartFlag()) return nullptr;
92  return *itStartCellHolder;
93  }
94 
95  private:
101  Weight getFinalCellState(ACellHolder* cellHolder,
102  const std::vector<WeightedRelation<ACellHolder>>& cellHolderRelations) const
103  {
104  AutomatonCell& cell = cellHolder->getAutomatonCell();
105 
106  // Throw if this cell has already been traversed in this recursion cycle
107  if (cell.hasCycleFlag()) {
108  B2DEBUG(100, "Cycle detected");
109  throw (CycleException());
110  }
111 
112  if (cell.hasAssignedFlag()) {
113  return cell.getCellState();
114 
115  } else {
116  // Mark cell in order to detect if it was already traversed in this recursion cycle
117  cell.setCycleFlag();
118 
119  Weight finalCellState = updateState(cellHolder, cellHolderRelations);
120 
121  // Unmark the cell
122  cell.unsetCycleFlag();
123  return finalCellState;
124  }
125  }
126 
128  Weight updateState(ACellHolder* cellHolder,
129  const std::vector<WeightedRelation<ACellHolder>>& cellHolderRelations) const
130  {
131  AutomatonCell& cell = cellHolder->getAutomatonCell();
132 
133  // --- blocked cells do not contribute a continuation ---
134  // Redundant check.
135  if (cell.hasMaskedFlag()) {
136  cell.setAssignedFlag();
137  return NAN;
138  }
139 
140  //--- Search for neighbors ---
141  Weight maxStateWithContinuation = NAN;
142 
143  // Flag to keep track whether the best continuation lies on a prioriy path
144  bool isPriorityPath = false;
145 
146  auto continuations = asRange(
147  std::equal_range(cellHolderRelations.begin(), cellHolderRelations.end(), cellHolder));
148 
149  // Check neighbors now
150  for (const WeightedRelation<ACellHolder>& relation : continuations) {
151  // Advance to valid neighbor
152  ACellHolder* neighborCellHolder = relation.getTo();
153 
154  // Skip dead ends (should not happen)
155  if (not neighborCellHolder) continue;
156 
157  AutomatonCell& neighborCell = neighborCellHolder->getAutomatonCell();
158  // Skip masked continuations
159  if (neighborCell.hasMaskedFlag()) continue;
160 
161  // Invalidate a possible start flag since the neighbor has an ancestors
162  neighborCell.unsetStartFlag();
163 
164  // Get the value of the neighbor
165  Weight neighborCellState = getFinalCellState(neighborCellHolder, cellHolderRelations);
166 
167  // Add the value of the connetion to the gain value
168  Weight stateWithContinuation = neighborCellState + relation.getWeight();
169 
170  // Remember only the maximum value of all neighbors
171  if (std::isnan(maxStateWithContinuation) or maxStateWithContinuation < stateWithContinuation) {
172  maxStateWithContinuation = stateWithContinuation;
173  // Remember whether the best continuation marks a priorty path
174  // construction ensures that priority paths have at least two elements
175  isPriorityPath = neighborCell.hasPriorityFlag() or neighborCell.hasPriorityPathFlag();
176  }
177  } // end for relations
178 
179  if (std::isnan(maxStateWithContinuation)) {
180  // No valid neighbor contributed a connection to the cell
181  maxStateWithContinuation = 0;
182  }
183 
184  // The value of this cell is only its own weight
185  maxStateWithContinuation += cell.getCellWeight();
186 
187  // Set the value
188  cell.setCellState(maxStateWithContinuation);
189  cell.setPriorityPathFlag(isPriorityPath);
190 
191  // Mark this cell as having its final value
192  cell.setAssignedFlag();
193 
194  // Return the just determined value
195  return cell.getCellState();
196  }
197 
198  private:
203  void prepareCellFlags(const std::vector<ACellHolder*>& cellHolders) const
204  {
205  for (ACellHolder* cellHolder : cellHolders) {
206  AutomatonCell& cell = cellHolder->getAutomatonCell();
207  cell.unsetTemporaryFlags();
208  if (cell.hasMaskedFlag()) continue;
209  cell.setCellState(NAN);
210  }
211  }
212  };
213  }
215 }
Belle2::TrackFindingCDC::AutomatonCell::hasAssignedFlag
bool hasAssignedFlag() const
Gets the current state of the already assigned marker flag.
Definition: AutomatonCell.h:150
Belle2::TrackFindingCDC::AutomatonCell::hasMaskedFlag
bool hasMaskedFlag() const
Gets the current state of the masked marker flag.
Definition: AutomatonCell.h:228
Belle2::TrackFindingCDC::CellularAutomaton::updateState
Weight updateState(ACellHolder *cellHolder, const std::vector< WeightedRelation< ACellHolder >> &cellHolderRelations) const
Updates the state of a cell considering all continuations recursively.
Definition: CellularAutomaton.h:136
Belle2::TrackFindingCDC::CellularAutomaton::getFinalCellState
Weight getFinalCellState(ACellHolder *cellHolder, const std::vector< WeightedRelation< ACellHolder >> &cellHolderRelations) const
Gets the cell state of the cell holder.
Definition: CellularAutomaton.h:109
Belle2::TrackFindingCDC::CellularAutomaton::prepareCellFlags
void prepareCellFlags(const std::vector< ACellHolder * > &cellHolders) const
Helper function to prepare the stats.
Definition: CellularAutomaton.h:211
Belle2::TrackFindingCDC::AutomatonCell::hasPriorityPathFlag
bool hasPriorityPathFlag() const
Gets the current state of the priority path marker flag.
Definition: AutomatonCell.h:186
Belle2::TrackFindingCDC::AutomatonCell::setStartFlag
void setStartFlag(bool setTo=true)
Sets the start marker flag to the given value. Default value true.
Definition: AutomatonCell.h:156
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFindingCDC::AutomatonCell
Cell used by the cellular automata.
Definition: AutomatonCell.h:39
Belle2::TrackFindingCDC::AutomatonCell::hasPriorityFlag
bool hasPriorityFlag() const
Gets the current state of the do not use flag marker flag.
Definition: AutomatonCell.h:305
Belle2::TrackFindingCDC::AutomatonCell::hasCycleFlag
bool hasCycleFlag() const
Gets the current state of the cycle marker flag.
Definition: AutomatonCell.h:204
Belle2::TrackFindingCDC::WeightedRelation
Type for two related objects with a weight.
Definition: CDCSegment2D.h:36
Belle2::TrackFindingCDC::CellularAutomaton::applyTo
ACellHolder * applyTo(const std::vector< ACellHolder * > &cellHolders, const std::vector< WeightedRelation< ACellHolder >> &cellHolderRelations) const
Applies the cellular automaton to the collection of cells and its neighborhood.
Definition: CellularAutomaton.h:53
Belle2::TrackFindingCDC::AutomatonCell::hasStartFlag
bool hasStartFlag() const
Gets the current state of the start marker flag.
Definition: AutomatonCell.h:168
Belle2::TrackFindingCDC::AutomatonCell::unsetStartFlag
void unsetStartFlag()
Resets the start marker flag to false.
Definition: AutomatonCell.h:162
Belle2::TrackFindingCDC::AutomatonCell::getCellState
Weight getCellState() const
Getter for the cell state.
Definition: AutomatonCell.h:106