Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-12-17 09:22:01

0001 /*!
0002  * \file Fun4AllSingleDstPileupInputManager.cc
0003  * \author Hugo Pereira Da Costa <hugo.pereira-da-costa@cea.fr>
0004  */
0005 
0006 #include "Fun4AllSingleDstPileupInputManager.h"
0007 #include "Fun4AllDstPileupMerger.h"
0008 
0009 #include <ffaobjects/RunHeader.h>
0010 
0011 #include <fun4all/DBInterface.h>
0012 #include <fun4all/Fun4AllReturnCodes.h>
0013 #include <fun4all/Fun4AllServer.h>
0014 
0015 #include <phool/PHCompositeNode.h>
0016 #include <phool/PHNodeIOManager.h>
0017 #include <phool/PHNodeIntegrate.h>
0018 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0019 #include <phool/PHRandomSeed.h>
0020 #include <phool/getClass.h>
0021 #include <phool/phool.h>  // for PHWHERE, PHReadOnly, PHRunTree
0022 
0023 #include <TSystem.h>
0024 
0025 #include <gsl/gsl_randist.h>
0026 
0027 #include <cassert>
0028 #include <iostream>  // for operator<<, basic_ostream, endl
0029 #include <utility>   // for pair
0030 
0031 //_____________________________________________________________________________
0032 Fun4AllSingleDstPileupInputManager::Fun4AllSingleDstPileupInputManager(const std::string &name, const std::string &nodename, const std::string &topnodename)
0033   : Fun4AllInputManager(name, nodename, topnodename)
0034 {
0035   // initialize random generator
0036   const uint seed = PHRandomSeed();
0037   m_rng.reset(gsl_rng_alloc(gsl_rng_mt19937));
0038   gsl_rng_set(m_rng.get(), seed);
0039 }
0040 
0041 //_____________________________________________________________________________
0042 int Fun4AllSingleDstPileupInputManager::fileopen(const std::string &filenam)
0043 {
0044   /*
0045   this is largely copied from fun4all/Fun4AllDstInputManager::fileopen
0046   with additional code to handle the background IManager
0047   */
0048 
0049   auto *se = Fun4AllServer::instance();
0050   if (IsOpen())
0051   {
0052     std::cout << "Closing currently open file "
0053               << FileName()
0054               << " and opening " << filenam << std::endl;
0055     fileclose();
0056   }
0057   FileName(filenam);
0058   m_fullfilename = DBInterface::instance()->location(FileName());
0059   if (Verbosity() > 0)
0060   {
0061     std::cout << Name() << ": opening file " << m_fullfilename << std::endl;
0062   }
0063   // sanity check - the IManager must be nullptr when this method is executed
0064   // if not something is very very wrong and we must not continue
0065   assert(!m_IManager);
0066 
0067   // first read the runnode if not disabled
0068   if (m_ReadRunTTree)
0069   {
0070     m_IManager.reset(new PHNodeIOManager(m_fullfilename, PHReadOnly, PHRunTree));
0071     if (m_IManager->isFunctional())
0072     {
0073       m_runNode = se->getNode(m_RunNode, TopNodeName());
0074       m_IManager->read(m_runNode);
0075 
0076       // get the current run number
0077       auto *runheader = findNode::getClass<RunHeader>(m_runNode, "RunHeader");
0078       if (runheader)
0079       {
0080         SetRunNumber(runheader->get_RunNumber());
0081       }
0082       // delete our internal copy of the runnode when opening subsequent files
0083       assert(!m_runNodeCopy);
0084       m_runNodeCopy.reset(new PHCompositeNode("RUNNODECOPY"));
0085       if (!m_runNodeSum)
0086       {
0087         m_runNodeSum.reset(new PHCompositeNode("RUNNODESUM"));
0088       }
0089 
0090       {
0091         // read run node using temporary node iomanager
0092         PHNodeIOManager(m_fullfilename, PHReadOnly, PHRunTree).read(m_runNodeCopy.get());
0093       }
0094 
0095       PHNodeIntegrate integrate;
0096       integrate.RunNode(m_runNode);
0097       integrate.RunSumNode(m_runNodeSum.get());
0098       // run recursively over internal run node copy and integrate objects
0099       PHNodeIterator mainIter(m_runNodeCopy.get());
0100       mainIter.forEach(integrate);
0101       // we do not need to keep the internal copy, keeping it would crate
0102       // problems in case a subsequent file does not contain all the
0103       // runwise objects from the previous file. Keeping this copy would then
0104       // integrate the missing object again with the old copy
0105       m_runNodeCopy.reset();
0106     }
0107     // DLW: move the delete outside the if block to cover the case where isFunctional() fails
0108     m_IManager.reset();
0109   }
0110 
0111   // create internal dst node
0112   if (!m_dstNodeInternal)
0113   {
0114     m_dstNodeInternal.reset(new PHCompositeNode("DST_INTERNAL"));
0115   }
0116 
0117   // update dst node from fun4all server
0118   m_dstNode = se->getNode(InputNode(), TopNodeName());
0119 
0120   // open file in both active and background input manager
0121   m_IManager.reset(new PHNodeIOManager(m_fullfilename, PHReadOnly));
0122   m_IManager_background.reset(new PHNodeIOManager(m_fullfilename, PHReadOnly));
0123   if (m_IManager->isFunctional())
0124   {
0125     IsOpen(1);
0126     m_ievent_thisfile = 0;
0127     setBranches();                // set branch selections
0128     AddToFileOpened(FileName());  // add file to the list of files which were opened
0129     return 0;
0130   }
0131 
0132   std::cout << PHWHERE << ": " << Name() << " Could not open file " << FileName() << std::endl;
0133   m_IManager.reset();
0134   m_IManager_background.reset();
0135   return -1;
0136 }
0137 
0138 //_____________________________________________________________________________
0139 int Fun4AllSingleDstPileupInputManager::run(const int nevents)
0140 {
0141   if (!IsOpen())
0142   {
0143     if (FileListEmpty())
0144     {
0145       if (Verbosity() > 0)
0146       {
0147         std::cout << Name() << ": No Input file open" << std::endl;
0148       }
0149       return -1;
0150     }
0151 
0152     if (OpenNextFile())
0153     {
0154       std::cout << Name() << ": No Input file from filelist opened" << std::endl;
0155       return -1;
0156     }
0157   }
0158 
0159   if (Verbosity() > 3)
0160   {
0161     std::cout << "Getting Event from " << Name() << std::endl;
0162   }
0163 
0164 readagain:
0165 
0166   // read main event to dstNode
0167   auto *dummy = m_IManager->read(m_dstNode);
0168   int ncount = 0;
0169   while (dummy)
0170   {
0171     ++ncount;
0172     if (nevents > 0 && ncount >= nevents)
0173     {
0174       break;
0175     }
0176     dummy = m_IManager->read(m_dstNode);
0177   }
0178 
0179   if (!dummy)
0180   {
0181     fileclose();
0182     if (!OpenNextFile())
0183     {
0184       // NOLINTNEXTLINE(hicpp-avoid-goto)
0185       goto readagain;
0186     }
0187     else
0188     {
0189       return -1;
0190     }
0191   }
0192 
0193   m_ievent_total += ncount;
0194   m_ievent_thisfile += ncount;
0195 
0196   // check events consistency
0197   if (m_ievent_thisfile != static_cast<int>(m_IManager->getEventNumber()))
0198   {
0199     std::cout << PHWHERE
0200               << " inconsistent event counters between inputmanager and ionode manager: "
0201               << " m_ievent_thisfile: " << m_ievent_thisfile
0202               << " m_IManager->getEventNumber(): " << m_IManager->getEventNumber()
0203               << std::endl;
0204     gSystem->Exit(1);
0205   }
0206 
0207   // check if the local SubsysReco discards this event
0208   if (RejectEvent() != Fun4AllReturnCodes::EVENT_OK)
0209   {
0210     std::cout << "Fun4AllSingleDstPileupInputManager::run - skipped event " << m_ievent_thisfile - 1 << std::endl;
0211     // NOLINTNEXTLINE(hicpp-avoid-goto)
0212     goto readagain;
0213   }
0214 
0215   // load relevant DST nodes to internal pointers
0216   if (Verbosity() > 0)
0217   {
0218     std::cout << "Fun4AllSingleDstPileupInputManager::run - loaded event " << m_ievent_thisfile - 1 << std::endl;
0219   }
0220 
0221   Fun4AllDstPileupMerger merger;
0222   merger.load_nodes(m_dstNode);
0223 
0224   // generate background collisions
0225   const double mu = m_collision_rate * m_time_between_crossings * 1e-9;
0226 
0227   const int min_crossing = m_tmin / m_time_between_crossings;
0228   const int max_crossing = m_tmax / m_time_between_crossings;
0229   int ievent_thisfile = m_ievent_thisfile;
0230   int neventsbackground = 0;
0231   for (int icrossing = min_crossing; icrossing <= max_crossing; ++icrossing)
0232   {
0233     const double crossing_time = m_time_between_crossings * icrossing;
0234     const int ncollisions = gsl_ran_poisson(m_rng.get(), mu);
0235     for (int icollision = 0; icollision < ncollisions; ++icollision)
0236     {
0237       // try read
0238       if (!m_IManager_background->read(m_dstNodeInternal.get(), ievent_thisfile))
0239       {
0240         break;
0241       }
0242 
0243       // merge
0244 
0245       if (Verbosity() > 0)
0246       {
0247         std::cout << "Fun4AllSingleDstPileupInputManager::run - merged background event " << ievent_thisfile << " time: " << crossing_time << std::endl;
0248       }
0249       merger.copy_background_event(m_dstNodeInternal.get(), crossing_time);
0250 
0251       ++neventsbackground;
0252       ++ievent_thisfile;
0253     }
0254   }
0255 
0256   // jump event counter to the last background accepted event
0257   if (neventsbackground > 0)
0258   {
0259     PushBackEvents(-neventsbackground);
0260   }
0261   return 0;
0262 }
0263 
0264 //_____________________________________________________________________________
0265 int Fun4AllSingleDstPileupInputManager::fileclose()
0266 {
0267   if (!IsOpen())
0268   {
0269     std::cout << Name() << ": fileclose: No Input file open" << std::endl;
0270     return -1;
0271   }
0272   m_IManager.reset();
0273   m_IManager_background.reset();
0274   IsOpen(0);
0275   UpdateFileList();
0276   return 0;
0277 }
0278 
0279 //_____________________________________________________________________________
0280 int Fun4AllSingleDstPileupInputManager::BranchSelect(const std::string &branch, const int iflag)
0281 {
0282   int myflag = iflag;
0283   // if iflag > 0 the branch is set to read
0284   // if iflag = 0, the branch is set to NOT read
0285   // if iflag < 0 the branchname is erased from our internal branch read map
0286   // this does not have any effect on phool yet
0287   if (myflag < 0)
0288   {
0289     std::map<const std::string, int>::iterator branchiter;
0290     branchiter = m_branchread.find(branch);
0291     if (branchiter != m_branchread.end())
0292     {
0293       m_branchread.erase(branchiter);
0294     }
0295     return 0;
0296   }
0297 
0298   if (myflag > 0)
0299   {
0300     if (Verbosity() > 1)
0301     {
0302       std::cout << "Setting Root Tree Branch: " << branch << " to read" << std::endl;
0303     }
0304     myflag = 1;
0305   }
0306   else
0307   {
0308     if (Verbosity() > 1)
0309     {
0310       std::cout << "Setting Root Tree Branch: " << branch << " to NOT read" << std::endl;
0311     }
0312   }
0313   m_branchread[branch] = myflag;
0314   return 0;
0315 }
0316 
0317 //_____________________________________________________________________________
0318 int Fun4AllSingleDstPileupInputManager::setBranches()
0319 {
0320   if (m_IManager)
0321   {
0322     if (!m_branchread.empty())
0323     {
0324       std::map<const std::string, int>::const_iterator branchiter;
0325       for (branchiter = m_branchread.begin(); branchiter != m_branchread.end(); ++branchiter)
0326       {
0327         m_IManager->selectObjectToRead(branchiter->first, branchiter->second);
0328         if (m_IManager_background)
0329         {
0330           m_IManager_background->selectObjectToRead(branchiter->first, branchiter->second);
0331         }
0332         if (Verbosity() > 0)
0333         {
0334           std::cout << branchiter->first << " set to " << branchiter->second << std::endl;
0335         }
0336       }
0337     }
0338   }
0339   else
0340   {
0341     std::cout << PHWHERE << " " << Name() << ": You can only call this function after a file has been opened" << std::endl;
0342     std::cout << "Do not worry, the branches will be set as soon as you open a file" << std::endl;
0343     return -1;
0344   }
0345   return 0;
0346 }
0347 
0348 //_____________________________________________________________________________
0349 void Fun4AllSingleDstPileupInputManager::Print(const std::string &what) const
0350 {
0351   if (what == "ALL" || what == "BRANCH")
0352   {
0353     // loop over the map and print out the content (name and location in memory)
0354     std::cout << "--------------------------------------" << std::endl
0355               << std::endl;
0356     std::cout << "List of selected branches in Fun4AllSingleDstPileupInputManager " << Name() << ":" << std::endl;
0357 
0358     std::map<const std::string, int>::const_iterator iter;
0359     for (iter = m_branchread.begin(); iter != m_branchread.end(); ++iter)
0360     {
0361       std::cout << iter->first << " is switched ";
0362       if (iter->second)
0363       {
0364         std::cout << "ON";
0365       }
0366       else
0367       {
0368         std::cout << "OFF";
0369       }
0370       std::cout << std::endl;
0371     }
0372   }
0373   if ((what == "ALL" || what == "PHOOL") && m_IManager)
0374   {
0375     // loop over the map and print out the content (name and location in memory)
0376     std::cout << "--------------------------------------" << std::endl
0377               << std::endl;
0378     std::cout << "PHNodeIOManager print in Fun4AllSingleDstPileupInputManager " << Name() << ":" << std::endl;
0379     m_IManager->print();
0380   }
0381   Fun4AllInputManager::Print(what);
0382   return;
0383 }
0384 
0385 //_____________________________________________________________________________
0386 int Fun4AllSingleDstPileupInputManager::PushBackEvents(const int i)
0387 {
0388   if (m_IManager)
0389   {
0390     unsigned EventOnDst = m_IManager->getEventNumber();
0391     EventOnDst -= static_cast<unsigned>(i);
0392     m_ievent_thisfile -= i;
0393     m_ievent_total -= i;
0394     m_IManager->setEventNumber(EventOnDst);
0395     return 0;
0396   }
0397   std::cout << PHWHERE << Name() << ": could not push back events, Imanager is NULL"
0398             << " probably the dst is not open yet (you need to call fileopen or run 1 event for lists)" << std::endl;
0399   return -1;
0400 }