12 #include <analysis/modules/ParticleListManipulator/ParticleListManipulatorModule.h>
16 #include <framework/datastore/DataStore.h>
17 #include <framework/datastore/StoreArray.h>
18 #include <framework/datastore/StoreObjPtr.h>
21 #include <framework/logging/Logger.h>
24 #include <analysis/dataobjects/Particle.h>
25 #include <analysis/dataobjects/ParticleList.h>
28 #include <analysis/DecayDescriptor/ParticleListName.h>
29 #include <analysis/utility/ValueIndexPairSorting.h>
53 setDescription(
"Manipulates ParticleLists: copies/merges/performs particle selection");
54 setPropertyFlags(c_ParallelProcessingCertified);
57 addParam(
"outputListName", m_outputListName,
"Output ParticleList name");
59 vector<string> defaultList;
60 addParam(
"inputListNames", m_inputListNames,
61 "list of input ParticleList names", defaultList);
63 addParam(
"cut", m_cutParameter,
"Selection criteria to be applied", std::string(
""));
65 addParam(
"variable", m_variableName,
"Variable which defines the best duplicate (see ``selectLowest`` for ordering)",
66 std::string(
"mdstIndex"));
67 addParam(
"preferLowest", m_preferLowest,
68 "If true, duplicate with lowest value of ``variable`` is accepted, otherwise higher one.",
true);
70 addParam(
"writeOut", m_writeOut,
71 "If true, the output ParticleList will be saved by RootOutput. If false, it will be ignored when writing the file.",
false);
75 m_isSelfConjugatedParticle =
false;
78 void ParticleListManipulatorModule::initialize()
83 bool valid = m_decaydescriptor.init(m_outputListName);
85 B2ERROR(
"ParticleListManipulatorModule::initialize Invalid output ParticleList name: " << m_outputListName);
90 m_pdgCode = mother->getPDGCode();
92 m_outputAntiListName = ParticleListName::antiParticleListName(m_outputListName);
93 m_isSelfConjugatedParticle = (m_outputListName == m_outputAntiListName);
96 for (
const std::string& listName : m_inputListNames) {
97 if (listName == m_outputListName) {
98 B2ERROR(
"ParticleListManipulatorModule: cannot copy Particles from " << listName <<
99 " to itself! Use applyCuts() (ParticleSelector module) instead.");
100 }
else if (!m_decaydescriptor.init(listName)) {
101 B2ERROR(
"Invalid input ParticleList name: " << listName);
109 particleList.registerInDataStore(flags);
110 if (!m_isSelfConjugatedParticle) {
112 antiParticleList.registerInDataStore(flags);
115 m_variable = Variable::Manager::Instance().getVariable(m_variableName);
117 B2ERROR(
"Variable '" << m_variableName <<
"' is not available in Variable::Manager!");
119 m_cut = Variable::Cut::compile(m_cutParameter);
122 void ParticleListManipulatorModule::event()
125 m_particlesInTheList.clear();
129 bool existingList = plist.
isValid();
134 plist->initialize(m_pdgCode, m_outputListName);
136 if (!m_isSelfConjugatedParticle) {
139 antiPlist->initialize(-1 * m_pdgCode, m_outputAntiListName);
141 antiPlist->bindAntiParticleList(*(plist));
147 for (
unsigned i = 0; i < plist->getListSize(); i++) {
148 const Particle* particle = plist->getParticle(i);
150 std::vector<int> idSeq;
151 fillUniqueIdentifier(particle, idSeq);
152 m_particlesInTheList.push_back(idSeq);
157 typedef std::pair<double, unsigned int> ValueIndexPair;
158 std::vector<ValueIndexPair> valueToIndex;
161 for (
const auto& inputListName : m_inputListNames) {
164 std::vector<int> fsParticles = inPList->getList(ParticleList::EParticleType::c_FlavorSpecificParticle,
false);
165 const std::vector<int>& scParticles = inPList->getList(ParticleList::EParticleType::c_SelfConjugatedParticle,
false);
166 const std::vector<int>& fsAntiParticles = inPList->getList(ParticleList::EParticleType::c_FlavorSpecificParticle,
true);
168 fsParticles.insert(fsParticles.end(), scParticles.begin(), scParticles.end());
169 fsParticles.insert(fsParticles.end(), fsAntiParticles.begin(), fsAntiParticles.end());
171 for (
int fsParticle : fsParticles) {
172 const Particle* part = particles[fsParticle];
174 if (m_cut->check(part)) {
175 valueToIndex.emplace_back(m_variable->function(part), part->getArrayIndex());
182 if (m_preferLowest) {
183 std::stable_sort(valueToIndex.begin(), valueToIndex.end(), ValueIndexPairSorting::lowerPair<ValueIndexPair>);
185 std::stable_sort(valueToIndex.begin(), valueToIndex.end(), ValueIndexPairSorting::higherPair<ValueIndexPair>);
189 for (
const auto& candidate : valueToIndex) {
190 const Particle* part = particles[candidate.second];
192 std::vector<int> idSeq;
193 fillUniqueIdentifier(part, idSeq);
194 bool uniqueSeq = isUnique(idSeq);
197 plist->addParticle(part);
198 m_particlesInTheList.push_back(idSeq);
203 void ParticleListManipulatorModule::fillUniqueIdentifier(
const Particle* p, std::vector<int>& idSequence)
205 idSequence.push_back(p->getPDGCode());
207 if (p->getNDaughters() == 0) {
208 idSequence.push_back(p->getMdstArrayIndex());
210 idSequence.push_back(p->getNDaughters());
211 auto daughters = p->getDaughters();
213 sort(daughters.begin(), daughters.end(), [](
const auto a,
const auto b) {
214 return a->getPDGCode() > b->getPDGCode();
217 for (
const auto& daughter : daughters)
218 fillUniqueIdentifier(daughter, idSequence);
222 bool ParticleListManipulatorModule::isUnique(
const std::vector<int>& idSeqOUT)
224 for (
const auto& idSeqIN : m_particlesInTheList) {
225 bool sameSeq = (idSeqIN == idSeqOUT);