Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:16:16

0001 #include "SingleMvtxPoolInput.h"
0002 #include "Fun4AllStreamingInputManager.h"
0003 #include "mvtx_pool.h"
0004 
0005 #include <ffarawobjects/MvtxFeeIdInfov1.h>
0006 #include <ffarawobjects/MvtxRawEvtHeaderv2.h>
0007 #include <ffarawobjects/MvtxRawHitContainerv1.h>
0008 #include <ffarawobjects/MvtxRawHitv1.h>
0009 #include <fun4all/Fun4AllUtils.h>
0010 #include "MvtxRawDefs.h"
0011 
0012 #include <frog/FROG.h>
0013 #include <phool/PHCompositeNode.h>
0014 #include <phool/PHNodeIterator.h>  // for PHNodeIterator
0015 #include <phool/getClass.h>
0016 #include <phool/phool.h>
0017 #include <phool/sphenix_constants.h>
0018 
0019 #include <Event/Event.h>
0020 #include <Event/EventTypes.h>
0021 #include <Event/Eventiterator.h>
0022 #include <Event/fileEventiterator.h>
0023 
0024 #include <cassert>
0025 #include <cmath>
0026 #include <memory>
0027 #include <set>
0028 
0029 SingleMvtxPoolInput::SingleMvtxPoolInput(const std::string &name)
0030   : SingleStreamingInput(name)
0031 {
0032   plist = new Packet *[2];
0033   m_rawHitContainerName = "MVTXRAWHIT";
0034 }
0035 
0036 SingleMvtxPoolInput::~SingleMvtxPoolInput()
0037 {
0038   delete[] plist;
0039   for (auto &iter : poolmap)
0040   {
0041     if (Verbosity() > 2)
0042     {
0043       std::cout << "deleting mvtx pool for id " << iter.first << std::endl;
0044     }
0045     delete (iter.second);
0046   }
0047 }
0048 
0049 void SingleMvtxPoolInput::FillPool(const uint64_t minBCO)
0050 {
0051   if (AllDone())  // no more files and all events read
0052   {
0053     return;
0054   }
0055   while (GetEventiterator() == nullptr)  // at startup this is a null pointer
0056   {
0057     if (!OpenNextFile())
0058     {
0059       AllDone(1);
0060       return;
0061     }
0062   }
0063 
0064   //  std::set<uint64_t> saved_beamclocks;
0065   while (GetSomeMoreEvents())
0066   {
0067     std::unique_ptr<Event> evt(GetEventiterator()->getNextEvent());
0068     while (!evt)
0069     {
0070       fileclose();
0071       if (!OpenNextFile())
0072       {
0073         AllDone(1);
0074         return;
0075       }
0076       evt.reset(GetEventiterator()->getNextEvent());
0077     }
0078     if (Verbosity() > 2)
0079     {
0080       std::cout << "Fetching next Event" << evt->getEvtSequence() << std::endl;
0081     }
0082     RunNumber(evt->getRunNumber());
0083     if (GetVerbosity() > 1)
0084     {
0085       evt->identify();
0086     }
0087     if (evt->getEvtType() != DATAEVENT)
0088     {
0089       m_NumSpecialEvents++;
0090       continue;
0091     }
0092     int EventSequence = evt->getEvtSequence();
0093     int npackets = evt->getPacketList(plist, 2);
0094 
0095     if (npackets > 2)
0096     {
0097       exit(1);
0098     }
0099     for (int i = 0; i < npackets; i++)
0100     {
0101       // Ignoring packet not from MVTX detector
0102       if (Verbosity() > 1)
0103       {
0104         plist[i]->identify();
0105       }
0106 
0107       if (poolmap.find(plist[i]->getIdentifier()) == poolmap.end())
0108       {
0109         if (Verbosity() > 1)
0110         {
0111           std::cout << "starting new mvtx pool for packet " << plist[i]->getIdentifier() << std::endl;
0112         }
0113         poolmap[plist[i]->getIdentifier()] = new mvtx_pool();
0114       }
0115       poolmap[plist[i]->getIdentifier()]->addPacket(plist[i]);
0116       delete plist[i];
0117     }
0118     for (auto &iter : poolmap)
0119     {
0120       mvtx_pool *pool = iter.second;
0121       int num_feeId = pool->get_feeidSet_size();
0122       if (Verbosity() > 1)
0123       {
0124         std::cout << "Number of feeid in RCDAQ events: " << num_feeId << " for packet "
0125                   << iter.first << std::endl;
0126       }
0127       if (num_feeId > 0)
0128       {
0129         for (int i_fee{0}; i_fee < num_feeId; ++i_fee)
0130         {
0131           // auto feeId = pool->iValue(i_fee, "FEEID");
0132           auto feeId = pool->get_feeid(i_fee);
0133           auto link = MvtxRawDefs::decode_feeid(feeId);
0134 
0135           //          auto hbfSize = plist[i]->iValue(feeId, "NR_HBF");
0136           // auto num_strobes_old = pool->iValue(feeId, "NR_STROBES");
0137           // auto num_L1Trgs_old = pool->iValue(feeId, "NR_PHYS_TRG");
0138           auto num_strobes = pool->get_strbSet_size(feeId);
0139           auto num_L1Trgs = pool->get_trgSet_size(feeId);
0140           for (int iL1 = 0; iL1 < num_L1Trgs; ++iL1)
0141           {
0142             //  auto l1Trg_bco = pool->lValue(feeId, iL1, "L1_IR_BCO");
0143             auto l1Trg_bco = pool->get_L1_IR_BCO(feeId, iL1);
0144             //            auto l1Trg_bc  = plist[i]->iValue(feeId, iL1, "L1_IR_BC");
0145             m_FeeGTML1BCOMap[feeId].insert(l1Trg_bco);
0146             gtmL1BcoSet.emplace(l1Trg_bco);
0147           }
0148           m_FeeStrobeMap[feeId] += num_strobes;
0149           for (int i_strb{0}; i_strb < num_strobes; ++i_strb)
0150           {
0151             // auto strb_detField = pool->iValue(feeId, i_strb, "TRG_DET_FIELD");
0152             // uint64_t strb_bco = pool->lValue(feeId, i_strb, "TRG_IR_BCO");
0153             // auto strb_bc = pool->iValue(feeId, i_strb, "TRG_IR_BC");
0154             // auto num_hits = pool->iValue(feeId, i_strb, "TRG_NR_HITS");
0155             auto strb_detField = pool->get_TRG_DET_FIELD(feeId, i_strb);
0156             uint64_t strb_bco = pool->get_TRG_IR_BCO(feeId, i_strb);
0157             auto strb_bc = pool->get_TRG_IR_BC(feeId, i_strb);
0158             auto num_hits = pool->get_TRG_NR_HITS(feeId, i_strb);
0159             m_BclkStack.insert(strb_bco);
0160             m_FEEBclkMap[feeId] = strb_bco;
0161 
0162             if (strb_bco < minBCO - m_NegativeBco)
0163             {
0164               continue;
0165             }
0166 
0167             if (Verbosity() > 4)
0168             {
0169               std::cout << "evtno: " << EventSequence << ", Fee: " << feeId;
0170               std::cout << " Layer: " << link.layer << " Stave: " << link.stave;
0171               std::cout << " GBT: " << link.gbtid << ", bco: 0x" << std::hex << strb_bco << std::dec;
0172               std::cout << ", n_hits: " << num_hits << std::endl;
0173             }
0174             auto hits = pool->get_hits(feeId, i_strb);
0175             for (auto &&hit : hits)
0176             {
0177               auto newhit = std::make_unique<MvtxRawHitv1>();
0178               newhit->set_bco(strb_bco);
0179               newhit->set_strobe_bc(strb_bc);
0180               newhit->set_chip_bc(hit->bunchcounter);
0181               newhit->set_layer_id(link.layer);
0182               newhit->set_stave_id(link.stave);
0183               newhit->set_chip_id(
0184                   MvtxRawDefs::gbtChipId_to_staveChipId[link.gbtid][hit->chip_id]);
0185               newhit->set_row(hit->row_pos);
0186               newhit->set_col(hit->col_pos);
0187               if (StreamingInputManager())
0188               {
0189                 StreamingInputManager()->AddMvtxRawHit(strb_bco, newhit.get());
0190               }
0191               m_MvtxRawHitMap[strb_bco].push_back(newhit.release());
0192             }
0193             if (StreamingInputManager())
0194             {
0195               StreamingInputManager()->AddMvtxFeeIdInfo(strb_bco, feeId, strb_detField);
0196             }
0197           }
0198         }
0199       }
0200     }
0201     // Assign L1 trg to Strobe windows data.
0202     for (auto &lv1Bco : gtmL1BcoSet)
0203     {
0204       auto it = m_BclkStack.lower_bound(lv1Bco);
0205       auto const strb_it = (it == m_BclkStack.begin()) ? (*it == lv1Bco ? it : m_BclkStack.cend()) : --it;
0206       if (strb_it != m_BclkStack.cend())
0207       {
0208         if (StreamingInputManager())
0209         {
0210           StreamingInputManager()->AddMvtxL1TrgBco(*strb_it, lv1Bco);
0211         }
0212       }
0213       else if (m_BclkStack.empty())
0214       {
0215         continue;
0216       }
0217       else
0218       {
0219         std::cout << "ERROR: lv1Bco: 0x" << std::hex << lv1Bco << std::dec
0220                   << " is less than minimun strobe bco 0x" << std::hex
0221                   << *m_BclkStack.begin() << std::dec << std::endl;
0222         // assert(0);
0223       }
0224     }
0225     gtmL1BcoSet.clear();
0226   }
0227 }
0228 
0229 void SingleMvtxPoolInput::Print(const std::string &what) const
0230 {
0231   // TODO: adapt to MVTX case
0232 
0233   if (what == "ALL" || what == "FEEBCLK")
0234   {
0235     for (auto bcliter : m_FEEBclkMap)
0236     {
0237       std::cout << "FEE" << bcliter.first << " bclk: 0x"
0238                 << std::hex << bcliter.second << std::dec << std::endl;
0239     }
0240   }
0241   if (what == "ALL" || what == "STORAGE")
0242   {
0243     for (const auto &bcliter : m_MvtxRawHitMap)
0244     {
0245       std::cout << "Beam clock 0x" << std::hex << bcliter.first << std::dec << std::endl;
0246       for (const auto &feeiter : bcliter.second)
0247       {
0248         std::cout << "fee: " << feeiter->get_stave_id()
0249                   << " at " << std::hex << feeiter << std::dec << std::endl;
0250       }
0251     }
0252   }
0253   if (what == "ALL" || what == "STACK")
0254   {
0255     for (const auto &iter : m_BclkStack)
0256     {
0257       std::cout << "stacked bclk: 0x" << std::hex << iter << std::dec << std::endl;
0258     }
0259   }
0260   if (what == "ALL" || what == "GET_NR_STROBES")
0261   {
0262     for (const auto &iter : m_FeeStrobeMap)
0263     {
0264       std::cout << "Total number of strobes for feeid: " << iter.first << ", " << iter.second << std::endl;
0265     }
0266   }
0267 }
0268 
0269 void SingleMvtxPoolInput::CleanupUsedPackets(const uint64_t bclk)
0270 {
0271   m_BclkStack.erase(m_BclkStack.begin(), m_BclkStack.upper_bound(bclk));
0272   for (auto it = m_MvtxRawHitMap.begin(); it != m_MvtxRawHitMap.end() && (it->first <= bclk); it = m_MvtxRawHitMap.erase(it))
0273   {
0274     for (const auto &rawhit : it->second)
0275     {
0276       delete rawhit;
0277     }
0278   }
0279   m_MvtxRawHitMap.erase(m_MvtxRawHitMap.begin(), m_MvtxRawHitMap.upper_bound(bclk));
0280   m_FeeStrobeMap.erase(m_FeeStrobeMap.begin(), m_FeeStrobeMap.upper_bound(bclk));
0281   for (auto &[feeid, gtmbcoset] : m_FeeGTML1BCOMap)
0282   {
0283     gtmbcoset.erase(gtmbcoset.begin(), gtmbcoset.upper_bound(bclk));
0284   }
0285 }
0286 
0287 bool SingleMvtxPoolInput::CheckPoolDepth(const uint64_t bclk)
0288 {
0289   // if (m_FEEBclkMap.size() < 10)
0290   // {
0291   //   std::cout << "not all FEEs in map: " << m_FEEBclkMap.size() << std::endl;
0292   //   return true;
0293   // }
0294   for (auto iter : m_FEEBclkMap)
0295   {
0296     if (Verbosity() > 2)
0297     {
0298       std::cout << iter.first << " my bclk 0x" << std::hex << iter.second
0299                 << " req: 0x" << bclk << std::dec << std::endl;
0300     }
0301     // equal case when we have more strobe with same bco
0302     // due not synchronization
0303     if (iter.second <= bclk)
0304     {
0305       if (Verbosity() > 1)
0306       {
0307         std::cout << "FEE " << iter.first << " beamclock 0x" << std::hex << iter.second
0308                   << " smaller than req bclk: 0x" << bclk << std::dec << std::endl;
0309       }
0310       return false;
0311     }
0312   }
0313   return true;
0314 }
0315 
0316 void SingleMvtxPoolInput::ClearCurrentEvent()
0317 {
0318   // called interactively, to get rid of the current event
0319   uint64_t currentbclk = *m_BclkStack.begin();
0320   //  std::cout << "clearing bclk 0x" << std::hex << currentbclk << std::dec << std::endl;
0321   CleanupUsedPackets(currentbclk);
0322   // m_BclkStack.erase(currentbclk);
0323   return;
0324 }
0325 
0326 bool SingleMvtxPoolInput::GetSomeMoreEvents()
0327 {
0328   if (AllDone())
0329   {
0330     return false;
0331   }
0332   if (m_MvtxRawHitMap.empty())
0333   {
0334     return true;
0335   }
0336   uint64_t lowest_bclk = m_MvtxRawHitMap.begin()->first;
0337   //  lowest_bclk += m_BcoRange;
0338   lowest_bclk += m_BcoRange;
0339   std::set<int> toerase;
0340   for (auto bcliter : m_FEEBclkMap)
0341   {
0342     if (bcliter.second <= lowest_bclk)
0343     {
0344       uint64_t highest_bclk = m_MvtxRawHitMap.rbegin()->first;
0345       if ((highest_bclk - m_MvtxRawHitMap.begin()->first) < MaxBclkDiff())
0346       {
0347         // std::cout << "FEE " << bcliter.first << " bclk: "
0348         //      << std::hex << bcliter.second << ", req: " << lowest_bclk
0349         //       << " low: 0x" <<  m_MvtxRawHitMap.begin()->first << ", high: " << highest_bclk << ", delta: " << std::dec << (highest_bclk-m_MvtxRawHitMap.begin()->first)
0350         //      << std::dec << std::endl;
0351         return true;
0352       }
0353       else
0354       {
0355         std::cout << PHWHERE << Name() << ": erasing FEE " << bcliter.first
0356                   << " with stuck bclk: " << std::hex << bcliter.second
0357                   << " current bco range: 0x" << m_MvtxRawHitMap.begin()->first
0358                   << ", to: 0x" << highest_bclk << ", delta: " << std::dec
0359                   << (highest_bclk - m_MvtxRawHitMap.begin()->first)
0360                   << std::dec << std::endl;
0361         toerase.insert(bcliter.first);
0362       }
0363     }
0364   }
0365   for (auto iter : toerase)
0366   {
0367     m_FEEBclkMap.erase(iter);
0368   }
0369   return false;
0370 
0371   //  if (CheckPoolDepth(m_MvtxRawHitMap.begin()->first))
0372   //  {
0373   //     if (m_MvtxRawHitMap.size() >= 200)
0374   //     {
0375   //       return false;
0376   //     }
0377   // //  }
0378   //   return true;
0379 }
0380 
0381 void SingleMvtxPoolInput::CreateDSTNode(PHCompositeNode *topNode)
0382 {
0383   PHNodeIterator iter(topNode);
0384   PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0385   if (!dstNode)
0386   {
0387     dstNode = new PHCompositeNode("DST");
0388     topNode->addNode(dstNode);
0389   }
0390   PHNodeIterator iterDst(dstNode);
0391   PHCompositeNode *detNode = dynamic_cast<PHCompositeNode *>(iterDst.findFirst("PHCompositeNode", "MVTX"));
0392   if (!detNode)
0393   {
0394     detNode = new PHCompositeNode("MVTX");
0395     dstNode->addNode(detNode);
0396   }
0397 
0398   MvtxRawEvtHeader *mvtxEH = findNode::getClass<MvtxRawEvtHeader>(detNode, m_rawEventHeaderName);
0399   if (!mvtxEH)
0400   {
0401     mvtxEH = new MvtxRawEvtHeaderv2();
0402     PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(mvtxEH, m_rawEventHeaderName, "PHObject");
0403     detNode->addNode(newNode);
0404   }
0405 
0406   MvtxRawHitContainer *mvtxhitcont = findNode::getClass<MvtxRawHitContainer>(detNode, m_rawHitContainerName);
0407   if (!mvtxhitcont)
0408   {
0409     mvtxhitcont = new MvtxRawHitContainerv1();
0410     PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(mvtxhitcont, m_rawHitContainerName, "PHObject");
0411     detNode->addNode(newNode);
0412   }
0413 }
0414 
0415 void SingleMvtxPoolInput::ConfigureStreamingInputManager()
0416 {
0417   auto [runnumber, segment] = Fun4AllUtils::GetRunSegment(*(GetFileList().begin()));
0418 
0419   if (m_readStrWidthFromDB)
0420   {
0421     m_strobeWidth = MvtxRawDefs::getStrobeLength(runnumber);
0422     if (std::isnan(m_strobeWidth))
0423     {
0424       std::cout << PHWHERE << "WARNING: Strobe length is not defined for run " << runnumber;
0425       std::cout << " neither in the OCDB or DAQ DB. Exiting SingleMvtxPoolInput." << std::endl;
0426       // std::cout << "Defaulting to 89 mus strobe length" << std::endl;
0427       // m_strobeWidth = 89.;
0428       exit(1);
0429     }
0430   }
0431 
0432   if (!m_mvtx_is_standalone)
0433   {
0434     if (m_strobeWidth > 88.)
0435     {
0436       m_BcoRange = 1000;
0437       m_NegativeBco = 1000;
0438     }
0439     else if (m_strobeWidth > 9 && m_strobeWidth < 11)
0440     {
0441       m_BcoRange = 500;
0442       m_NegativeBco = 500;
0443     }
0444     else if (m_strobeWidth < 1)  // triggered mode
0445     {
0446       m_BcoRange = 3;
0447       m_NegativeBco = 0;
0448       if (StreamingInputManager())
0449       {
0450         StreamingInputManager()->runMvtxTriggered(true);
0451       }
0452     }
0453     else  // catchall for anyting else to set to a range based on the rhic clock
0454     {
0455       m_BcoRange = std::ceil(m_strobeWidth * 1000. / sphenix_constants::time_between_crossings);
0456       m_NegativeBco = std::ceil(m_strobeWidth * 1000. / sphenix_constants::time_between_crossings);
0457     }
0458   }
0459 
0460   if (Verbosity() > 1)
0461   {
0462     std::cout << "Mvtx strobe length " << m_strobeWidth << std::endl;
0463     std::cout << "Mvtx BCO range and negative bco range set based on strobe length " << m_BcoRange << ", " << m_NegativeBco << std::endl;
0464   }
0465 
0466   if (StreamingInputManager())
0467   {
0468     StreamingInputManager()->SetMvtxBcoRange(m_BcoRange);
0469     StreamingInputManager()->SetMvtxNegativeBco(m_NegativeBco);
0470   }
0471   return;
0472 }