Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:18:11

0001 #include "PHG4TruthInfoContainer.h"
0002 
0003 #include "PHG4Particle.h"
0004 #include "PHG4Shower.h"
0005 #include "PHG4VtxPoint.h"
0006 
0007 #include <boost/tuple/tuple.hpp>
0008 
0009 #include <limits>
0010 #include <string>
0011 
0012 using namespace std;
0013 
0014 PHG4TruthInfoContainer::PHG4TruthInfoContainer()
0015   : particlemap()
0016   , vtxmap()
0017   , particle_embed_flags()
0018   , vertex_embed_flags()
0019 {
0020 }
0021 
0022 PHG4TruthInfoContainer::~PHG4TruthInfoContainer() { Reset(); }
0023 
0024 void PHG4TruthInfoContainer::Reset()
0025 {
0026   for (auto& iter : particlemap)
0027   {
0028     delete iter.second;
0029   }
0030   particlemap.clear();
0031 
0032   for (auto& iter : sPHENIXprimaryparticlemap)
0033   {
0034     delete iter.second;
0035   }
0036   sPHENIXprimaryparticlemap.clear();
0037 
0038   for (auto& iter : vtxmap)
0039   {
0040     delete iter.second;
0041   }
0042   vtxmap.clear();
0043 
0044   for (auto& iter : showermap)
0045   {
0046     delete iter.second;
0047   }
0048   showermap.clear();
0049 
0050   particle_embed_flags.clear();
0051   vertex_embed_flags.clear();
0052 
0053   return;
0054 }
0055 
0056 void PHG4TruthInfoContainer::identify(ostream& os) const
0057 {
0058   os << "---particlemap--------------------------" << endl;
0059   for (auto iter : particlemap)
0060   {
0061     os << "particle id " << iter.first << endl;
0062     (iter.second)->identify();
0063   }
0064 
0065   os << "---vtxmap-------------------------------" << endl;
0066   for (auto vter : vtxmap)
0067   {
0068     os << "vtx id: " << vter.first << endl;
0069     (vter.second)->identify();
0070   }
0071 
0072   os << "---showermap-------------------------------" << endl;
0073   for (auto ster : showermap)
0074   {
0075     os << "shower id: " << ster.first << endl;
0076     (ster.second)->identify();
0077   }
0078 
0079   os << "---list of embeded track flags-------------------" << endl;
0080   for (auto particle_embed_flag : particle_embed_flags)
0081   {
0082     os << "embeded track id: " << particle_embed_flag.first
0083        << " flag: " << particle_embed_flag.second << endl;
0084   }
0085 
0086   os << "---list of embeded vtx flags-------------------" << endl;
0087   for (auto vertex_embed_flag : vertex_embed_flags)
0088   {
0089     os << "embeded vertex id: " << vertex_embed_flag.first
0090        << " flag: " << vertex_embed_flag.second << endl;
0091   }
0092 
0093   os << "---primary vertex-------------------" << endl;
0094   os << "Vertex " << GetPrimaryVertexIndex() << " is identified as the primary vertex" << endl;
0095 
0096   return;
0097 }
0098 
0099 PHG4TruthInfoContainer::ConstIterator
0100 PHG4TruthInfoContainer::AddParticle(const int trackid, PHG4Particle* newparticle)
0101 {
0102   int key = trackid;
0103   ConstIterator it;
0104   bool added = false;
0105   boost::tie(it, added) = particlemap.insert(std::make_pair(key, newparticle));
0106   if (added)
0107   {
0108     return it;
0109   }
0110 
0111   cerr << "PHG4TruthInfoContainer::AddParticle"
0112        << " - Attempt to add particle with existing trackid "
0113        << trackid << ": " << newparticle->get_name() << " id "
0114        << newparticle->get_track_id()
0115        << ", p = [" << newparticle->get_px()
0116        << ", " << newparticle->get_py()
0117        << ", " << newparticle->get_pz() << "], "
0118        << " parent ID " << newparticle->get_parent_id()
0119        << std::endl;
0120   return particlemap.end();
0121 }
0122 
0123 PHG4TruthInfoContainer::ConstIterator
0124 PHG4TruthInfoContainer::AddsPHENIXPrimaryParticle(const int trackid, PHG4Particle* newparticle)
0125 {
0126   int key = trackid;
0127   ConstIterator it;
0128   bool added = false;
0129 
0130   boost::tie(it, added) = sPHENIXprimaryparticlemap.insert(std::make_pair(key, newparticle));
0131   if (added)
0132   {
0133     return it;
0134   }
0135 
0136   cerr << "PHG4TruthInfoContainer::AddsPHENIXPrimaryParticle"
0137        << " - Attempt to add sPHENIX primary particle with existing trackid "
0138        << trackid << ": " << newparticle->get_name() << " id "
0139        << newparticle->get_track_id()
0140        << ", p = [" << newparticle->get_px()
0141        << ", " << newparticle->get_py()
0142        << ", " << newparticle->get_pz() << "], "
0143        << " parent ID " << newparticle->get_parent_id()
0144        << std::endl;
0145   return sPHENIXprimaryparticlemap.end();
0146 }
0147 
0148 PHG4Particle* PHG4TruthInfoContainer::GetParticle(const int trackid)
0149 {
0150   int key = trackid;
0151   Iterator it = particlemap.find(key);
0152   if (it != particlemap.end())
0153   {
0154     return it->second;
0155   }
0156   return nullptr;
0157 }
0158 
0159 PHG4Particle* PHG4TruthInfoContainer::GetParticle(const int trackid) const
0160 {
0161   int key = trackid;
0162   ConstIterator it = particlemap.find(key);
0163   if (it != particlemap.end())
0164   {
0165     return it->second;
0166   }
0167   return nullptr;
0168 }
0169 
0170 PHG4Particle* PHG4TruthInfoContainer::GetPrimaryParticle(const int trackid)
0171 {
0172   if (trackid <= 0)
0173   {
0174     return nullptr;
0175   }
0176   Iterator it = particlemap.find(trackid);
0177   if (it != particlemap.end())
0178   {
0179     return it->second;
0180   }
0181   return nullptr;
0182 }
0183 
0184 PHG4Particle* PHG4TruthInfoContainer::GetsPHENIXPrimaryParticle(const int trackid)
0185 {
0186   Iterator it = sPHENIXprimaryparticlemap.find(trackid);
0187   if (it != sPHENIXprimaryparticlemap.end())
0188   {
0189     return it->second;
0190   }
0191   return nullptr;
0192 }
0193 
0194 PHG4VtxPoint* PHG4TruthInfoContainer::GetVtx(const int vtxid)
0195 {
0196   int key = vtxid;
0197   VtxIterator it = vtxmap.find(key);
0198   if (it != vtxmap.end())
0199   {
0200     return it->second;
0201   }
0202   return nullptr;
0203 }
0204 
0205 PHG4VtxPoint* PHG4TruthInfoContainer::GetPrimaryVtx(const int vtxid)
0206 {
0207   if (vtxid <= 0)
0208   {
0209     return nullptr;
0210   }
0211   VtxIterator it = vtxmap.find(vtxid);
0212   if (it != vtxmap.end())
0213   {
0214     return it->second;
0215   }
0216   return nullptr;
0217 }
0218 
0219 PHG4Shower* PHG4TruthInfoContainer::GetShower(const int showerid)
0220 {
0221   int key = showerid;
0222   ShowerIterator it = showermap.find(key);
0223   if (it != showermap.end())
0224   {
0225     return it->second;
0226   }
0227   return nullptr;
0228 }
0229 
0230 PHG4Shower* PHG4TruthInfoContainer::GetPrimaryShower(const int showerid)
0231 {
0232   if (showerid <= 0)
0233   {
0234     return nullptr;
0235   }
0236   ShowerIterator it = showermap.find(showerid);
0237   if (it != showermap.end())
0238   {
0239     return it->second;
0240   }
0241   return nullptr;
0242 }
0243 
0244 PHG4TruthInfoContainer::ConstVtxIterator
0245 PHG4TruthInfoContainer::AddVertex(const int id, PHG4VtxPoint* newvtx)
0246 {
0247   int key = id;
0248   ConstVtxIterator it;
0249   bool added = false;
0250 
0251   if (vtxmap.find(id) != vtxmap.end())
0252   {
0253     cout << "trying to add existing vtx " << id
0254          << " vtx pos: " << endl;
0255     (vtxmap.find(id)->second)->identify();
0256     identify();
0257   }
0258 
0259   boost::tie(it, added) = vtxmap.insert(std::make_pair(key, newvtx));
0260   if (added)
0261   {
0262     newvtx->set_id(key);
0263     return it;
0264   }
0265 
0266   cerr << "PHG4TruthInfoContainer::AddVertex"
0267        << " - Attempt to add vertex with existing id " << id << std::endl;
0268   return vtxmap.end();
0269 }
0270 
0271 PHG4TruthInfoContainer::ConstShowerIterator
0272 PHG4TruthInfoContainer::AddShower(const int id, PHG4Shower* newshower)
0273 {
0274   int key = id;
0275   ConstShowerIterator it;
0276   bool added = false;
0277 
0278   if (showermap.find(id) != showermap.end())
0279   {
0280     cout << "trying to add existing shower " << id
0281          << " shower pos: " << endl;
0282     (showermap.find(id)->second)->identify();
0283     identify();
0284   }
0285 
0286   boost::tie(it, added) = showermap.insert(std::make_pair(key, newshower));
0287   if (added)
0288   {
0289     newshower->set_id(key);
0290     return it;
0291   }
0292 
0293   cerr << "PHG4TruthInfoContainer::AddShower"
0294        << " - Attempt to add shower with existing id " << id << std::endl;
0295   return showermap.end();
0296 }
0297 
0298 int PHG4TruthInfoContainer::maxtrkindex() const
0299 {
0300   int key = 0;
0301   if (!particlemap.empty())
0302   {
0303     key = particlemap.rbegin()->first;
0304   }
0305   if (key < 0)
0306   {
0307     key = 0;
0308   }
0309   return key;
0310 }
0311 
0312 int PHG4TruthInfoContainer::mintrkindex() const
0313 {
0314   int key = 0;
0315   if (!particlemap.empty())
0316   {
0317     key = particlemap.begin()->first;
0318   }
0319   if (key > 0)
0320   {
0321     key = 0;
0322   }
0323   return key;
0324 }
0325 
0326 int PHG4TruthInfoContainer::maxvtxindex() const
0327 {
0328   int key = 0;
0329   if (!vtxmap.empty())
0330   {
0331     key = vtxmap.rbegin()->first;
0332   }
0333   if (key < 0)
0334   {
0335     key = 0;
0336   }
0337   return key;
0338 }
0339 
0340 int PHG4TruthInfoContainer::minvtxindex() const
0341 {
0342   int key = 0;
0343   if (!vtxmap.empty())
0344   {
0345     key = vtxmap.begin()->first;
0346   }
0347   if (key > 0)
0348   {
0349     key = 0;
0350   }
0351   return key;
0352 }
0353 
0354 int PHG4TruthInfoContainer::maxshowerindex() const
0355 {
0356   int key = 0;
0357   if (!showermap.empty())
0358   {
0359     key = showermap.rbegin()->first;
0360   }
0361   if (key < 0)
0362   {
0363     key = 0;
0364   }
0365   return key;
0366 }
0367 
0368 int PHG4TruthInfoContainer::minshowerindex() const
0369 {
0370   int key = 0;
0371   if (!showermap.empty())
0372   {
0373     key = showermap.begin()->first;
0374   }
0375   if (key > 0)
0376   {
0377     key = 0;
0378   }
0379   return key;
0380 }
0381 
0382 void PHG4TruthInfoContainer::delete_particle(Iterator piter)
0383 {
0384   delete piter->second;
0385   particlemap.erase(piter);
0386   return;
0387 }
0388 
0389 void PHG4TruthInfoContainer::delete_particle(int trackid)
0390 {
0391   Iterator it = particlemap.find(trackid);
0392   if (it != particlemap.end())
0393   {
0394     delete_particle(it);
0395   }
0396 }
0397 
0398 void PHG4TruthInfoContainer::delete_vtx(VtxIterator viter)
0399 {
0400   delete viter->second;
0401   vtxmap.erase(viter);
0402   return;
0403 }
0404 
0405 void PHG4TruthInfoContainer::delete_vtx(int vtxid)
0406 {
0407   VtxIterator it = vtxmap.find(vtxid);
0408   if (it != vtxmap.end())
0409   {
0410     delete_vtx(it);
0411   }
0412 }
0413 
0414 void PHG4TruthInfoContainer::delete_shower(ShowerIterator siter)
0415 {
0416   delete siter->second;
0417   showermap.erase(siter);
0418   return;
0419 }
0420 
0421 int PHG4TruthInfoContainer::isEmbeded(const int trackid) const
0422 {
0423   //I think here for G4 secondarys we just check their primary's embedding ID
0424   int trackid_embed = trackid;
0425   if (trackid_embed <= 0)
0426   {
0427     const PHG4Particle* p = GetParticle(trackid_embed);
0428     if(p) trackid_embed = p->get_primary_id();
0429   }
0430   std::map<int, int>::const_iterator iter = particle_embed_flags.find(trackid_embed);
0431   if (iter == particle_embed_flags.end())
0432   {
0433     return 0;
0434   }
0435 
0436   return iter->second;
0437 }
0438 
0439 int PHG4TruthInfoContainer::isEmbededVtx(const int vtxid) const
0440 {
0441   std::map<int, int>::const_iterator iter = vertex_embed_flags.find(vtxid);
0442   if (iter == vertex_embed_flags.end())
0443   {
0444     return 0;
0445   }
0446 
0447   return iter->second;
0448 }
0449 
0450 bool PHG4TruthInfoContainer::is_primary_vtx(const PHG4VtxPoint* v) const
0451 {
0452   return (v->get_id() > 0);
0453 }
0454 
0455 bool PHG4TruthInfoContainer::is_primary(const PHG4Particle* p) const
0456 {
0457   return (p->get_track_id() > 0);
0458 }
0459 
0460 // this is O(log N)...
0461 bool PHG4TruthInfoContainer::is_sPHENIX_primary(const PHG4Particle* p) const
0462 {
0463   // sPHENIX primary particles are in the sPHENIXprimaryparticlemap
0464   return (sPHENIXprimaryparticlemap.find(p->get_track_id()) != sPHENIXprimaryparticlemap.end());
0465 }
0466 
0467 int PHG4TruthInfoContainer::GetPrimaryVertexIndex() const
0468 {
0469   ConstVtxRange vrange = GetPrimaryVtxRange();
0470 
0471   int highest_embedding_ID = numeric_limits<int>::min();
0472   int vtx_id_for_highest_embedding_ID = 0;
0473   for (auto iter = vrange.first; iter != vrange.second; ++iter)
0474   {
0475     //        cout <<"PHG4TruthInfoContainer::GetPrimaryVertexIndex - vertex ID "<<iter->first<<" embedding ID "<< isEmbededVtx(iter->first) <<": ";
0476     // iter->second->identify();
0477     const int embedding_ID = isEmbededVtx(iter->first);
0478 
0479     if (embedding_ID >= highest_embedding_ID)
0480     {
0481       highest_embedding_ID = embedding_ID;
0482       vtx_id_for_highest_embedding_ID = iter->first;
0483     }
0484   }
0485 
0486   if (highest_embedding_ID == numeric_limits<int>::min())
0487   {
0488     cout << "PHG4TruthInfoContainer::GetPrimaryVertexIndex - "
0489          << "WARNING: no valid primary vertex. Return an invalid ID of 0"
0490          << endl;
0491     return 0;
0492   }
0493 
0494   return vtx_id_for_highest_embedding_ID;
0495 }
0496 
0497 bool operator==(const PHG4TruthInfoContainer::Map::value_type& lhs, const PHG4TruthInfoContainer::Map::value_type& rhs)
0498 {
0499   return *lhs.second == *rhs.second;
0500 }
0501 
0502 bool operator==(const PHG4TruthInfoContainer::VtxMap::value_type& lhs, const PHG4TruthInfoContainer::VtxMap::value_type& rhs)
0503 {
0504   return *lhs.second == *rhs.second;
0505 }
0506 
0507 bool operator==(const PHG4TruthInfoContainer::ShowerMap::value_type& lhs, const PHG4TruthInfoContainer::ShowerMap::value_type& rhs)
0508 {
0509   return *lhs.second == *rhs.second;
0510 }