packages/mvtx/MvtxCombinedRawDataDecoder.cc

0001 /*!
0002  * \file MvtxCombinedRawDataDecoder.cc
0003  * \author Cameron Dean <cameron.dean@cern.ch>
0004  * \author Jakub Kvapil <jakub.kvapil@cern.ch>
0005  */
0006
0007 #include "MvtxCombinedRawDataDecoder.h"
0008
0009 #include <fun4allraw/MvtxRawDefs.h>
0010 #include <trackbase/MvtxDefs.h>
0011 #include <trackbase/MvtxEventInfov3.h>
0012 #include <trackbase/TrkrHitSet.h>
0013 #include <trackbase/TrkrHitSetContMvtxHelperv1.h>
0014 #include <trackbase/TrkrHitSetContainerv1.h>
0015 #include <trackbase/TrkrHitv2.h>
0016
0017 #include <fun4all/Fun4AllServer.h>
0018
0019 #include <ffarawobjects/Gl1Packet.h>
0020 #include <ffarawobjects/Gl1RawHit.h>
0021 #include <ffarawobjects/MvtxFeeIdInfov1.h>
0022 #include <ffarawobjects/MvtxRawEvtHeaderv1.h>
0023 #include <ffarawobjects/MvtxRawEvtHeaderv2.h>
0024 #include <ffarawobjects/MvtxRawHitContainerv1.h>
0025 #include <ffarawobjects/MvtxRawHitv1.h>
0026
0027 #include <fun4all/Fun4AllReturnCodes.h>
0028
0029 #include <phool/PHCompositeNode.h>
0030 #include <phool/PHNodeIterator.h>
0031 #include <phool/getClass.h>
0032 #include <phool/recoConsts.h>
0033 #include <phool/sphenix_constants.h>
0034
0035 #include <cdbobjects/CDBTTree.h>
0036 #include <ffamodules/CDBInterface.h>  // for accessing the MVTX hot pixel file from the CDB
0037
0038 #include <cassert>
0039 #include <iterator>
0040
0041 namespace
0042 {
0043   std::string MvtxHitSetHelperName("TRKR_MVTXHITSETHELPER");
0044 }
0045
0046 //_________________________________________________________
0047 MvtxCombinedRawDataDecoder::MvtxCombinedRawDataDecoder(const std::string &name)
0048   : SubsysReco(name)
0049 {
0050 }
0051
0052 //___________________________________________________________________________
0053 void MvtxCombinedRawDataDecoder::CreateNodes(PHCompositeNode *topNode)
0054 {
0055   // get dst node
0056   PHNodeIterator iter(topNode);
0057   PHCompositeNode *dstNode =
0058       dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
0059   if (!dstNode)
0060   {
0061     std::cout << "MvtxCombinedRawDataDecoder::InitRun - DST Node missing, "
0062                  "doing nothing."
0063               << std::endl;
0064     exit(1);
0065   }
0066
0067   // create mvtx hitset helper container if needed
0068   mvtx_hit_set_helper =
0069       findNode::getClass<TrkrHitSetContMvtxHelper>(topNode, MvtxHitSetHelperName);
0070   if (!mvtx_hit_set_helper)
0071   {
0072     // find or create TRKR node
0073     auto *trkrNode = dynamic_cast<PHCompositeNode *>(
0074         iter.findFirst("PHCompositeNode", "TRKR"));
0075     if (!trkrNode)
0076     {
0077       trkrNode = new PHCompositeNode("TRKR");
0078       dstNode->addNode(trkrNode);
0079     }
0080
0081     // create container and add to the tree
0082     mvtx_hit_set_helper = new TrkrHitSetContMvtxHelperv1;
0083     auto *newNode = new PHIODataNode<PHObject>(mvtx_hit_set_helper, MvtxHitSetHelperName,
0084                                                "PHObject");
0085     trkrNode->addNode(newNode);
0086   }
0087
0088   // create hitset container if needed
0089   hit_set_container =
0090       findNode::getClass<TrkrHitSetContainer>(topNode, "TRKR_HITSET");
0091   if (!hit_set_container)
0092   {
0093     // find or create TRKR node
0094     auto *trkrNode = dynamic_cast<PHCompositeNode *>(
0095         iter.findFirst("PHCompositeNode", "TRKR"));
0096     if (!trkrNode)
0097     {
0098       trkrNode = new PHCompositeNode("TRKR");
0099       dstNode->addNode(trkrNode);
0100     }
0101
0102     // create container and add to the tree
0103     hit_set_container = new TrkrHitSetContainerv1;
0104     auto *newNode = new PHIODataNode<PHObject>(hit_set_container, "TRKR_HITSET",
0105                                                "PHObject");
0106     trkrNode->addNode(newNode);
0107   }
0108
0109   // Check if MVTX event header already exists
0110   auto *mvtxNode = dynamic_cast<PHCompositeNode *>(
0111       iter.findFirst("PHCompositeNode", "MVTX"));
0112   if (!mvtxNode)
0113   {
0114     mvtxNode = new PHCompositeNode("MVTX");
0115     dstNode->addNode(mvtxNode);
0116   }
0117
0118   mvtx_event_header =
0119       findNode::getClass<MvtxEventInfo>(mvtxNode, "MVTXEVENTHEADER");
0120   if (!mvtx_event_header)
0121   {
0122     mvtx_event_header = new MvtxEventInfov3();
0123     auto *newHeader = new PHIODataNode<PHObject>(
0124         mvtx_event_header, "MVTXEVENTHEADER", "PHObject");
0125     mvtxNode->addNode(newHeader);
0126   }
0127
0128   mvtx_raw_event_header =
0129       findNode::getClass<MvtxRawEvtHeader>(topNode, m_MvtxRawEvtHeaderNodeName);
0130
0131   mvtx_raw_hit_container =
0132       findNode::getClass<MvtxRawHitContainer>(topNode, m_MvtxRawHitNodeName);
0133 }
0134
0135 //_____________________________________________________________________
0136 void MvtxCombinedRawDataDecoder::GetNodes(PHCompositeNode *topNode)
0137 {
0138   mvtx_raw_event_header =
0139       findNode::getClass<MvtxRawEvtHeader>(topNode, m_MvtxRawEvtHeaderNodeName);
0140   if (Verbosity() >= 3)
0141   {
0142     mvtx_raw_event_header->identify();
0143   }
0144
0145   mvtx_raw_hit_container =
0146       findNode::getClass<MvtxRawHitContainer>(topNode, m_MvtxRawHitNodeName);
0147   if (Verbosity() >= 3)
0148   {
0149     mvtx_raw_hit_container->identify();
0150   }
0151
0152   hit_set_container =
0153       findNode::getClass<TrkrHitSetContainer>(topNode, "TRKR_HITSET");
0154
0155   mvtx_hit_set_helper =
0156       findNode::getClass<TrkrHitSetContMvtxHelper>(topNode, MvtxHitSetHelperName);
0157
0158   mvtx_event_header =
0159       findNode::getClass<MvtxEventInfo>(topNode, "MVTXEVENTHEADER");
0160
0161   // Could we just get the first strobe BCO instead of setting this to 0?
0162   // Possible problem, what if the first BCO isn't the mean, then we'll shift tracker hit sets? Probably not a bad thing but depends on hit stripping
0163   //  uint64_t gl1rawhitbco = gl1 ? gl1->get_bco() : 0;
0164   auto *gl1 = findNode::getClass<Gl1Packet>(topNode, "GL1RAWHIT");
0165   if (gl1)
0166   {
0167     gl1rawhitbco = gl1->lValue(0, "BCO");
0168   }
0169   else
0170   {
0171     auto *oldgl1 = findNode::getClass<Gl1RawHit>(topNode, "GL1RAWHIT");
0172     if (oldgl1)
0173     {
0174       gl1rawhitbco = oldgl1->get_bco();
0175     }
0176   }
0177   if (gl1rawhitbco == 0 && (Verbosity() >= 4))
0178   {
0179     std::cout << PHWHERE << "Could not get gl1 raw hit" << std::endl;
0180   }
0181 }
0182
0183 //_____________________________________________________________________
0184 int MvtxCombinedRawDataDecoder::Init(PHCompositeNode * /*topNode*/)
0185 {
0186   return Fun4AllReturnCodes::EVENT_OK;
0187 }
0188
0189 //____________________________________________________________________________..
0190 int MvtxCombinedRawDataDecoder::InitRun(PHCompositeNode *topNode)
0191 {
0192   CreateNodes(topNode);
0193
0194   if (!mvtx_raw_event_header || !mvtx_raw_hit_container)
0195   {
0196     Fun4AllServer::instance()->unregisterSubsystem(this);
0197     std::cout << PHWHERE << "::" << __func__ << ": Could not get \""
0198               << m_MvtxRawHitNodeName << " or " << m_MvtxRawEvtHeaderNodeName << "\" from Node Tree" << std::endl;
0199     std::cout << "Unregistering subsystem and continuing on" << std::endl;
0200
0201   }
0202   if (dynamic_cast<MvtxRawEvtHeaderv1 *>(mvtx_raw_event_header))
0203   {
0204     std::cout << PHWHERE <<  "MvtxCombinedRawDataDecoder::GetNodes() !!!WARNING!!! using obsolete MvtxRawEvtHeaderv1.";
0205     std::cout << " Unregistering MvtxCombinedRawDataDecoder SubsysReco." << std::endl;
0206     Fun4AllServer::instance()->unregisterSubsystem(this);
0207   }
0208
0209   auto *rc = recoConsts::instance();
0210   int runNumber = rc->get_IntFlag("RUNNUMBER");
0211
0212   if (m_readStrWidthFromDB)
0213   {
0214     m_strobeWidth = MvtxRawDefs::getStrobeLength(runNumber);
0215   }
0216   if (std::isnan(m_strobeWidth))
0217   {
0218     std::cout << "MvtxCombinedRawDataDecoder::InitRun - strobe width is undefined for this run, defaulting to 89 mus" << std::endl;
0219     m_strobeWidth = 89;
0220   }
0221   if (m_strobeWidth < 1)
0222   {
0223     runMvtxTriggered(true);
0224   }
0225   //! save MVTX strobe width
0226   rc->set_FloatFlag("MvtxStrobeWidth", m_strobeWidth);
0227
0228   // Load the hot pixel map from the CDB
0229   if (m_doOfflineMasking)
0230   {
0231     m_hot_pixel_mask = new MvtxPixelMask();
0232     m_hot_pixel_mask->load_from_CDB();
0233   }
0234
0235   return Fun4AllReturnCodes::EVENT_OK;
0236 }
0237
0238 //___________________________________________________________________________
0239 int MvtxCombinedRawDataDecoder::process_event(PHCompositeNode *topNode)
0240 {
0241   GetNodes(topNode);
0242
0243   // get the last 40 bits by bit shifting left then right to match
0244   // to the mvtx bco
0245   auto lbshift = gl1rawhitbco << 24U;
0246   auto gl1bco = lbshift >> 24U;
0247
0248   // std::vector<std::pair<uint64_t, uint32_t> > strobe_bc_pairs;
0249   // std::set<uint64_t> l1BCOs = mvtx_raw_event_header->getMvtxLvL1BCO();
0250   // auto mvtxbco = *l1BCOs.begin();
0251   // if (gl1rawhitbco && (Verbosity() > 2))
0252   // {
0253   //   std::cout << "MVTX header BCO " << mvtxbco << " and GL1 BCO " << gl1bco
0254   //             << std::endl;
0255   // }
0256
0257   for (const auto &L1 : mvtx_raw_event_header->getMvtxLvL1BCO())
0258   {
0259     mvtx_event_header->add_L1_BCO(L1);
0260   }
0261   auto nMvtxFeeIdInfo = mvtx_raw_event_header->get_nFeeIdInfo();
0262   for (size_t i{}; i < nMvtxFeeIdInfo; ++i)
0263   {
0264     auto *feeIdInfo = mvtx_raw_event_header->get_feeIdInfo(i);
0265     mvtx_event_header->add_strobe_BCO(feeIdInfo->get_bco());
0266   }
0267   if (Verbosity() >= 3)
0268   {
0269     mvtx_raw_hit_container->identify();
0270   }
0271   const auto &strobe_list = mvtx_event_header->get_strobe_BCOs();
0272   auto it_strb_bco_zero = strobe_list.upper_bound(gl1bco);
0273   auto str_wGL1_idx = std::distance(strobe_list.cbegin(), it_strb_bco_zero) - 1;
0274
0275   uint64_t hit_strobe = -1;  // Initialise to -1 for debugging
0276   uint8_t layer = 0;
0277   uint8_t stave = 0;
0278   uint8_t chip = 0;
0279   uint16_t row = 0;
0280   uint16_t col = 0;
0281
0282   for (unsigned int i = 0; i < mvtx_raw_hit_container->get_nhits(); i++)
0283   {
0284     mvtx_rawhit = mvtx_raw_hit_container->get_hit(i);
0285     hit_strobe = mvtx_rawhit->get_bco();
0286     layer = mvtx_rawhit->get_layer_id();
0287     stave = mvtx_rawhit->get_stave_id();
0288     chip = mvtx_rawhit->get_chip_id();
0289     row = mvtx_rawhit->get_row();
0290     col = mvtx_rawhit->get_col();
0291
0292     int strobe_index = -20;
0293     const auto it = strobe_list.find(hit_strobe);
0294     if (it == strobe_list.cend())
0295     {
0296       std::cout << "Warning: hit strobe BCO " << hit_strobe << " is not found in evet combined strobe list:" << std::endl;
0297       for (const auto &strobe : strobe_list)
0298       {
0299         std::cout << "0x" << std::hex << strobe << std::dec << std::endl;
0300       }
0301     }
0302     else
0303     {
0304       strobe_index = static_cast<int>(std::distance(strobe_list.cbegin(), it) - str_wGL1_idx);
0305     }
0306
0307     // int bcodiff = gl1rawhitbco ? strobe - gl1bco : 0;
0308     //   double timeElapsed = bcodiff * 0.1065;  // 106 ns rhic clock
0309     //   int index = m_mvtx_is_triggered ? 0 : std::ceil(timeElapsed / m_strobeWidth);
0310
0311     if (strobe_index < -16 || strobe_index > 15)
0312     {
0313       std::cout << "Strobe index: " << strobe_index << " out of range" << std::endl;
0314       continue;  // Index is out of the 5-bit signed range
0315     }
0316
0317     if (Verbosity() >= 10)
0318     {
0319       mvtx_rawhit->identify();
0320     }
0321
0322     const TrkrDefs::hitsetkey hitsetkey =
0323         MvtxDefs::genHitSetKey(layer, stave, chip, strobe_index);
0324     if (!hitsetkey)
0325     {
0326       continue;
0327     }
0328
0329     mvtx_hit_set_helper->addHitSetKey(strobe_index, hitsetkey);
0330
0331     // get matching hitset
0332     const auto hitset_it = hit_set_container->findOrAddHitSet(hitsetkey);
0333
0334     // generate hit key
0335     const TrkrDefs::hitkey hitkey = MvtxDefs::genHitKey(col, row);
0336
0337     // find existing hit, or create
0338     auto *hit = hitset_it->second->getHit(hitkey);
0339     if (hit)
0340     {
0341       if (Verbosity() > 1)
0342       {
0343         std::cout << PHWHERE << "::" << __func__
0344                   << " - duplicated hit, hitsetkey: " << hitsetkey
0345                   << " hitkey: " << hitkey << std::endl;
0346       }
0347       continue;
0348     }
0349
0350     if (m_doOfflineMasking)
0351     {
0352       if (!m_hot_pixel_mask->is_masked(mvtx_rawhit))
0353       {  // Check if the pixel is masked
0354         hit = new TrkrHitv2;
0355         hitset_it->second->addHitSpecificKey(hitkey, hit);
0356       }
0357     }
0358     else
0359     {
0360       hit = new TrkrHitv2;
0361       hitset_it->second->addHitSpecificKey(hitkey, hit);
0362     }
0363   }
0364
0365   return Fun4AllReturnCodes::EVENT_OK;
0366 }
0367
0368 //_____________________________________________________________________
0369 int MvtxCombinedRawDataDecoder::End(PHCompositeNode * /*topNode*/)
0370 {
0371   return Fun4AllReturnCodes::EVENT_OK;
0372 }
0373
0374
0375 // void MvtxCombinedRawDataDecoder::removeDuplicates(
0376 //     std::vector<std::pair<uint64_t, uint32_t> > &v)
0377 // {
0378 //   auto end = v.end();
0379 //   for (auto it = v.begin(); it != end; ++it)
0380 //   {
0381 //     end = remove(it + 1, end, *it);
0382 //   }
0383 //   v.erase(end, v.end());
0384 // }