Belle II Software  light-2205-abys
ParticleListManipulatorModule.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 // Own include
10 #include <analysis/modules/ParticleListManipulator/ParticleListManipulatorModule.h>
11 
12 
13 // framework - DataStore
14 #include <framework/datastore/DataStore.h>
15 
16 // framework aux
17 #include <framework/logging/Logger.h>
18 
19 // utilities
20 #include <analysis/DecayDescriptor/ParticleListName.h>
21 #include <analysis/utility/ValueIndexPairSorting.h>
22 
23 using namespace std;
24 using namespace Belle2;
25 
26 //-----------------------------------------------------------------
27 // Register module
28 //-----------------------------------------------------------------
29 
30 REG_MODULE(ParticleListManipulator);
31 
32 //-----------------------------------------------------------------
33 // Implementation
34 //-----------------------------------------------------------------
35 
36 ParticleListManipulatorModule::ParticleListManipulatorModule():
37  m_variable(nullptr)
38 
39 {
40  setDescription("Manipulates ParticleLists: copies/merges/performs particle selection");
42 
43  // Add parameters
44  addParam("outputListName", m_outputListName, "Output ParticleList name");
45 
46  vector<string> defaultList;
47  addParam("inputListNames", m_inputListNames,
48  "list of input ParticleList names", defaultList);
49 
50  addParam("cut", m_cutParameter, "Selection criteria to be applied", std::string(""));
51 
52  addParam("variable", m_variableName, "Variable which defines the best duplicate (see ``selectLowest`` for ordering)",
53  std::string("mdstIndex"));
54  addParam("preferLowest", m_preferLowest,
55  "If true, duplicate with lowest value of ``variable`` is accepted, otherwise higher one.", true);
56 
57  addParam("writeOut", m_writeOut,
58  "If true, the output ParticleList will be saved by RootOutput. If false, it will be ignored when writing the file.", false);
59 
60  // initializing the rest of private members
61  m_pdgCode = 0;
63 }
64 
66 {
67  m_particles.isRequired();
68 
69  m_pdgCode = 0;
70 
71  // check the validity of output ParticleList name
73  if (!valid)
74  B2ERROR("ParticleListManipulatorModule::initialize Invalid output ParticleList name: " << m_outputListName);
75 
76  // Output particle
78 
79  m_pdgCode = mother->getPDGCode();
80 
81  // Some labels are reserved for the particle loader which loads all particles of the corresponding type.
82  // If people apply cuts or merge particle lists, the resulting particle lists are not allowed to have these names.
83  // Otherwise, very dangerous bugs could be introduced.
84  string listLabel = mother->getLabel();
85  // For final state particles we protect the label "all".
86  if (Const::finalStateParticlesSet.contains(Const::ParticleType(abs(m_pdgCode))) and listLabel == "all") {
87  B2FATAL("You have tried to create the list " << m_outputListName <<
88  " but the label 'all' is forbidden for user-defined lists of final-state particles." <<
89  " It could introduce *very* dangerous bugs.");
90  } else if ((listLabel == "MC") or (listLabel == "V0" and not(("K_S0:mdst" == m_inputListNames[0])
91  or ("Lambda0:mdst" == m_inputListNames[0]) or ("gamma:v0mdst" == m_inputListNames[0])))) {
92  // the labels MC, and V0 are also protected
93  // copying of some B2BII V0 lists has to be allowed to not break the FEI
94  B2FATAL("You have tried to create the list " << m_outputListName <<
95  " but the label " << listLabel << " is not allowed for merged or copied particle lists.");
96  }
97 
100 
101  // Input lists
102  for (const std::string& listName : m_inputListNames) {
103  if (listName == m_outputListName) {
104  B2ERROR("ParticleListManipulatorModule: cannot copy Particles from " << listName <<
105  " to itself! Use applyCuts() (ParticleSelector module) instead.");
106  } else if (!m_decaydescriptor.init(listName)) {
107  B2ERROR("Invalid input ParticleList name: " << listName);
108  } else {
110  }
111  }
112 
114  m_particleList.registerInDataStore(m_outputListName, flags);
116  m_antiParticleList.registerInDataStore(m_outputAntiListName, flags);
117  }
118 
120  if (!m_variable) {
121  B2ERROR("Variable '" << m_variableName << "' is not available in Variable::Manager!");
122  }
123  if (!(m_variable->variabletype == Variable::Manager::VariableDataType::c_double
124  or m_variable->variabletype == Variable::Manager::VariableDataType::c_int)) {
125  B2ERROR("Variable '" << m_variableName << "' has wrong data type! It must be either double or integer.");
126  }
128 }
129 
131 {
132  // clear the list
133  m_particlesInTheList.clear();
134 
135  bool existingList = m_particleList.isValid();
136 
137  if (!existingList) {
138  // new particle list: create it
139  m_particleList.create();
141 
143  m_antiParticleList.create();
145 
146  m_antiParticleList->bindAntiParticleList(*(m_particleList));
147  }
148  } else {
149  // output list already contains Particles
150  // fill m_particlesInTheList with unique
151  // identifiers of particles already in
152  for (unsigned i = 0; i < m_particleList->getListSize(); i++) {
153  const Particle* particle = m_particleList->getParticle(i);
154 
155  std::vector<int> idSeq;
156  fillUniqueIdentifier(particle, idSeq);
157  m_particlesInTheList.push_back(idSeq);
158  }
159  }
160 
161  // create list of candidate indices and corresponding sorting values
162  typedef std::pair<double, unsigned int> ValueIndexPair;
163  std::vector<ValueIndexPair> valueToIndex;
164 
165  // fill all particles from input lists that pass selection criteria into comparison list
166  for (const auto& inputListName : m_inputListNames) {
167  const StoreObjPtr<ParticleList> inPList(inputListName);
168  if (!inPList.isValid()) continue;
169 
170  std::vector<int> fsParticles = inPList->getList(ParticleList::EParticleType::c_FlavorSpecificParticle, false);
171  const std::vector<int>& scParticles = inPList->getList(ParticleList::EParticleType::c_SelfConjugatedParticle, false);
172  const std::vector<int>& fsAntiParticles = inPList->getList(ParticleList::EParticleType::c_FlavorSpecificParticle, true);
173 
174  fsParticles.insert(fsParticles.end(), scParticles.begin(), scParticles.end());
175  fsParticles.insert(fsParticles.end(), fsAntiParticles.begin(), fsAntiParticles.end());
176 
177  for (int fsParticle : fsParticles) {
178  const Particle* part = m_particles[fsParticle];
179 
180  if (m_cut->check(part)) {
181  double value = std::numeric_limits<double>::quiet_NaN();;
182  if (std::holds_alternative<double>(m_variable->function(part))) {
183  value = std::get<double>(m_variable->function(part));
184  } else if (std::holds_alternative<int>(m_variable->function(part))) {
185  value = std::get<int>(m_variable->function(part));
186  }
187  valueToIndex.emplace_back(value, part->getArrayIndex());
188  }
189  }
190  }
191 
192  // use stable sort to make sure we keep the relative order of elements with
193  // same value as it was before
194  if (m_preferLowest) {
195  std::stable_sort(valueToIndex.begin(), valueToIndex.end(), ValueIndexPairSorting::lowerPair<ValueIndexPair>);
196  } else {
197  std::stable_sort(valueToIndex.begin(), valueToIndex.end(), ValueIndexPairSorting::higherPair<ValueIndexPair>);
198  }
199 
200  // starting from the best candidate add all particles to output list that are not already in it
201  for (const auto& candidate : valueToIndex) {
202  const Particle* part = m_particles[candidate.second];
203 
204  std::vector<int> idSeq;
205  fillUniqueIdentifier(part, idSeq);
206  bool uniqueSeq = isUnique(idSeq);
207 
208  if (uniqueSeq) {
209  m_particleList->addParticle(part);
210  m_particlesInTheList.push_back(idSeq);
211  }
212  }
213 }
214 
215 void ParticleListManipulatorModule::fillUniqueIdentifier(const Particle* p, std::vector<int>& idSequence)
216 {
217  idSequence.push_back(p->getPDGCode());
218 
219  if (p->getNDaughters() == 0) {
220  idSequence.push_back(p->getMdstArrayIndex());
221  } else {
222  idSequence.push_back(p->getNDaughters());
223  auto daughters = p->getDaughters();
224  // sorting the daughters by their pdgCode to identify decay chains only differing by the order of their daughters
225  sort(daughters.begin(), daughters.end(), [](const auto a, const auto b) {
226  return a->getPDGCode() > b->getPDGCode();
227  });
228  // this is not FSP (go one level down)
229  for (const auto& daughter : daughters)
230  fillUniqueIdentifier(daughter, idSequence);
231  }
232 }
233 
234 bool ParticleListManipulatorModule::isUnique(const std::vector<int>& idSeqOUT)
235 {
236  for (const auto& idSeqIN : m_particlesInTheList) {
237  bool sameSeq = (idSeqIN == idSeqOUT);
238  if (sameSeq)
239  return false;
240  }
241 
242  return true;
243 }
The ParticleType class for identifying different particle types.
Definition: Const.h:289
static const ParticleSet finalStateParticlesSet
set of final set particles that can be created by the ParticleLoader
Definition: Const.h:538
EStoreFlags
Flags describing behaviours of objects etc.
Definition: DataStore.h:69
@ c_WriteOut
Object/array should be saved by output modules.
Definition: DataStore.h:70
@ c_DontWriteOut
Object/array should be NOT saved by output modules.
Definition: DataStore.h:71
Represents a particle in the DecayDescriptor.
int getPDGCode() const
Return PDG code.
std::string getLabel() const
The label of this particle, "default" returned, when no label set.
bool init(const std::string &str)
Initialise the DecayDescriptor from given string.
const DecayDescriptorParticle * getMother() const
return mother.
static std::unique_ptr< GeneralCut > compile(const std::string &cut)
Creates an instance of a cut and returns a unique_ptr to it, if you need a copy-able object instead y...
Definition: GeneralCut.h:84
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:208
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition: Module.h:80
bool m_isSelfConjugatedParticle
flag that indicates whether an anti-particle does not exist and therefore the output anti-ParticleLis...
bool m_preferLowest
Select the duplicate with the lowest value (instead of highest).
virtual void initialize() override
Initialize the Module.
virtual void event() override
Event processor.
std::string m_variableName
Variable which defines the best duplicate.
std::vector< std::string > m_inputListNames
input ParticleList names
StoreArray< Particle > m_particles
StoreArray of Particles.
std::vector< std::vector< int > > m_particlesInTheList
This vector holds unique identifiers (vector of ints) of all particles that are already included in t...
std::unique_ptr< Variable::Cut > m_cut
cut object which performs the cuts
bool isUnique(const std::vector< int > &idSequence)
Compares input idSequence to all sequences already included in the list.
const Variable::Manager::Var * m_variable
Variable which defines the best duplicate selection.
DecayDescriptor m_decaydescriptor
Decay descriptor of the particle being selected.
StoreObjPtr< ParticleList > m_particleList
output particle list
std::string m_outputAntiListName
output anti-particle list name
StoreObjPtr< ParticleList > m_antiParticleList
output anti-particle list
void fillUniqueIdentifier(const Particle *p, std::vector< int > &idSequence)
Fills unique identifier for the input particle.
std::string m_outputListName
output ParticleList name
Class to store reconstructed particles.
Definition: Particle.h:74
int getArrayIndex() const
Returns this object's array index (in StoreArray), or -1 if not found.
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:95
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:110
const Var * getVariable(std::string name)
Get the variable belonging to the given key.
Definition: Manager.cc:57
static Manager & Instance()
get singleton instance.
Definition: Manager.cc:25
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:560
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:650
std::string antiParticleListName(const std::string &listName)
Returns name of anti-particle-list corresponding to listName.
Abstract base class for different kinds of events.
Definition: ClusterUtils.h:23
VariableDataType variabletype
data type of variable
Definition: Manager.h:133
FunctionPtr function
Pointer to function.
Definition: Manager.h:147