Back to home page

sPhenix code displayed by LXR
 
 

    

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

0001 //____________________________________________________________________________..
0002 //
0003 // This is a template for a Fun4All SubsysReco module with all methods from the
0004 // $OFFLINE_MAIN/include/fun4all/SubsysReco.h baseclass
0005 // You do not have to implement all of them, you can just remove unused methods
0006 // here and in dNdEtaINTT.h.
0007 //
0008 // dNdEtaINTT(const std::string &name = "dNdEtaINTT")
0009 // everything is keyed to dNdEtaINTT, duplicate names do work but it makes
0010 // e.g. finding culprits in logs difficult or getting a pointer to the module
0011 // from the command line
0012 //
0013 // dNdEtaINTT::~dNdEtaINTT()
0014 // this is called when the Fun4AllServer is deleted at the end of running. Be
0015 // mindful what you delete - you do loose ownership of object you put on the
0016 // node tree
0017 //
0018 // int dNdEtaINTT::Init(PHCompositeNode *topNode)
0019 // This method is called when the module is registered with the Fun4AllServer.
0020 // You can create historgrams here or put objects on the node tree but be aware
0021 // that modules which haven't been registered yet did not put antyhing on the
0022 // node tree
0023 //
0024 // int dNdEtaINTT::InitRun(PHCompositeNode *topNode)
0025 // This method is called when the first event is read (or generated). At
0026 // this point the run number is known (which is mainly interesting for raw data
0027 // processing). Also all objects are on the node tree in case your module's
0028 // action depends on what else is around. Last chance to put nodes under the DST
0029 // Node We mix events during readback if branches are added after the first
0030 // event
0031 //
0032 // int dNdEtaINTT::process_event(PHCompositeNode *topNode)
0033 // called for every event. Return codes trigger actions, you find them in
0034 // $OFFLINE_MAIN/include/fun4all/Fun4AllReturnCodes.h
0035 //   everything is good:
0036 //     return Fun4AllReturnCodes::EVENT_OK
0037 //   abort event reconstruction, clear everything and process next event:
0038 //     return Fun4AllReturnCodes::ABORT_EVENT;
0039 //   proceed but do not save this event in output (needs output manager
0040 //   setting):
0041 //     return Fun4AllReturnCodes::DISCARD_EVENT;
0042 //   abort processing:
0043 //     return Fun4AllReturnCodes::ABORT_RUN
0044 // all other integers will lead to an error and abort of processing
0045 //
0046 // int dNdEtaINTT::ResetEvent(PHCompositeNode *topNode)
0047 // If you have internal data structures (arrays, stl containers) which needs
0048 // clearing after each event, this is the place to do that. The nodes under the
0049 // DST node are cleared by the framework
0050 //
0051 // int dNdEtaINTT::EndRun(const int runnumber)
0052 // This method is called at the end of a run when an event from a new run is
0053 // encountered. Useful when analyzing multiple runs (raw data). Also called at
0054 // the end of processing (before the End() method)
0055 //
0056 // int dNdEtaINTT::End(PHCompositeNode *topNode)
0057 // This is called at the end of processing. It needs to be called by the macro
0058 // by Fun4AllServer::End(), so do not forget this in your macro
0059 //
0060 // int dNdEtaINTT::Reset(PHCompositeNode *topNode)
0061 // not really used - it is called before the dtor is called
0062 //
0063 // void dNdEtaINTT::Print(const std::string &what) const
0064 // Called from the command line - useful to print information when you need it
0065 //
0066 //____________________________________________________________________________..
0067 
0068 #include "dNdEtaINTT.h"
0069 
0070 #include <fun4all/Fun4AllReturnCodes.h>
0071 #include <fun4all/PHTFileServer.h>
0072 #include <intt/CylinderGeomIntt.h>
0073 #include <phool/PHCompositeNode.h>
0074 #include <trackbase/TrkrHitv2.h>
0075 
0076 #include <limits>
0077 #include <sstream>
0078 
0079 namespace
0080 {
0081 template <class T> void CleanVec(std::vector<T> &v)
0082 {
0083     std::vector<T>().swap(v);
0084     v.shrink_to_fit();
0085 }
0086 } // namespace
0087 
0088 //____________________________________________________________________________..
0089 dNdEtaINTT::dNdEtaINTT(const std::string &name, const std::string &outputfile, const bool &isData)
0090     : SubsysReco(name)
0091     , _get_hepmc_info(true)
0092     , _get_truth_cluster(true)
0093     , _get_reco_cluster(true)
0094     , _get_centrality(true)
0095     , _get_intt_data(true)
0096     , _get_inttrawhit(false)
0097     , _get_trkr_hit(true)
0098     , _get_phg4_info(true)
0099     , _get_pmt_info(true)
0100     , _get_trigger_info(true)
0101     , _outputFile(outputfile)
0102     , IsData(isData)
0103     , eventheader(nullptr)
0104     , m_geneventmap(nullptr)
0105     , m_genevt(nullptr)
0106     , svtx_evalstack(nullptr)
0107     , truth_eval(nullptr)
0108     , clustereval(nullptr)
0109     , hiteval(nullptr)
0110     , dst_clustermap(nullptr)
0111     , clusterhitmap(nullptr)
0112     , inttrawhitcontainer(nullptr)
0113     , hitsets(nullptr)
0114     , _tgeometry(nullptr)
0115     , _intt_geom_container(nullptr)
0116     , m_truth_info(nullptr)
0117     , m_CentInfo(nullptr)
0118 {
0119     if (Verbosity() >= VERBOSITY_MORE)
0120         std::cout << "dNdEtaINTT::dNdEtaINTT(const std::string &name) Calling ctor" << std::endl;
0121 }
0122 
0123 //____________________________________________________________________________..
0124 dNdEtaINTT::~dNdEtaINTT()
0125 {
0126     if (Verbosity() >= VERBOSITY_MORE)
0127         std::cout << "dNdEtaINTT::~dNdEtaINTT() Calling dtor" << std::endl;
0128 }
0129 
0130 //____________________________________________________________________________..
0131 int dNdEtaINTT::Init(PHCompositeNode *topNode)
0132 {
0133     std::cout << "dNdEtaINTT::Init(PHCompositeNode *topNode) Initializing" << std::endl << "Running on Data or simulation? -> IsData = " << IsData << std::endl;
0134 
0135     PHTFileServer::get().open(_outputFile, "RECREATE");
0136     outtree = new TTree("EventTree", "EventTree");
0137     outtree->Branch("event", &event_);
0138     // outtree->Branch("event_counter", &evt_sequence_);
0139     if (_get_centrality)
0140     {
0141         if (!IsData)
0142         {
0143             outtree->Branch("ncoll", &ncoll_);
0144             outtree->Branch("npart", &npart_);
0145             outtree->Branch("centrality_bimp", &centrality_bimp_);
0146             outtree->Branch("centrality_impactparam", &centrality_impactparam_);
0147         }
0148         // outtree->Branch("event", &event_);
0149         outtree->Branch("clk", &clk);
0150         outtree->Branch("femclk", &femclk);
0151         outtree->Branch("is_min_bias", &is_min_bias);
0152         outtree->Branch("is_min_bias_wozdc", &is_min_bias_wozdc);
0153         outtree->Branch("MBD_centrality", &centrality_mbd_);
0154         outtree->Branch("MBD_z_vtx", &mbd_z_vtx);
0155         outtree->Branch("MBD_south_npmt", &mbd_south_npmt);
0156         outtree->Branch("MBD_north_npmt", &mbd_north_npmt);
0157         outtree->Branch("MBD_south_charge_sum", &mbd_south_charge_sum);
0158         outtree->Branch("MBD_north_charge_sum", &mbd_north_charge_sum);
0159         outtree->Branch("MBD_charge_sum", &mbd_charge_sum);
0160         outtree->Branch("MBD_charge_asymm", &mbd_charge_asymm);
0161         outtree->Branch("MBD_nhitsoverths_south", &mbd_nhitsoverths_south);
0162         outtree->Branch("MBD_nhitsoverths_north", &mbd_nhitsoverths_north);
0163         outtree->Branch("MBD_PMT_charge", &m_pmt_q, "MBD_PMT_charge[128]/F");
0164         // if (_get_pmt_info)
0165         // {
0166         //     for (unsigned int i = 0; i < 128; i++)
0167         //     {
0168         //         std::ostringstream pmtNumber;
0169         //         pmtNumber << std::setfill('0') << std::setw(3) << std::to_string(i);
0170         //         std::string branchName = "MBD_PMT" + pmtNumber.str() + "_charge";
0171         //         outtree->Branch(branchName.c_str(), &m_pmt_q[i]);
0172         //     }
0173         // }
0174     }
0175     if (_get_intt_data)
0176     {
0177         // outtree->Branch("event_counter", &event_);
0178         outtree->Branch("INTT_BCO", &intt_bco);
0179 
0180         if (!IsData)
0181         {
0182             outtree->Branch("NTruthVtx", &NTruthVtx_);
0183             outtree->Branch("TruthPV_trig_x", &TruthPV_trig_x_);
0184             outtree->Branch("TruthPV_trig_y", &TruthPV_trig_y_);
0185             outtree->Branch("TruthPV_trig_z", &TruthPV_trig_z_);
0186             // HepMC::GenParticle informaiton (final states, after boost and rotation)
0187             outtree->Branch("NHepMCFSPart", &NHepMCFSPart_);
0188             outtree->Branch("signal_process_id", &signal_process_id_);
0189             outtree->Branch("HepMCFSPrtl_Pt", &HepMCFSPrtl_Pt_);
0190             outtree->Branch("HepMCFSPrtl_Eta", &HepMCFSPrtl_Eta_);
0191             outtree->Branch("HepMCFSPrtl_Phi", &HepMCFSPrtl_Phi_);
0192             outtree->Branch("HepMCFSPrtl_E", &HepMCFSPrtl_E_);
0193             outtree->Branch("HepMCFSPrtl_PID", &HepMCFSPrtl_PID_);
0194             outtree->Branch("HepMCFSPrtl_prodx", &HepMCFSPrtl_prodx_);
0195             outtree->Branch("HepMCFSPrtl_prody", &HepMCFSPrtl_prody_);
0196             outtree->Branch("HepMCFSPrtl_prodz", &HepMCFSPrtl_prodz_);
0197             // Primary PHG4Particle information
0198             outtree->Branch("NPrimaryG4P", &NPrimaryG4P_);
0199             outtree->Branch("PrimaryG4P_Pt", &PrimaryG4P_Pt_);
0200             outtree->Branch("PrimaryG4P_Eta", &PrimaryG4P_Eta_);
0201             outtree->Branch("PrimaryG4P_Phi", &PrimaryG4P_Phi_);
0202             outtree->Branch("PrimaryG4P_E", &PrimaryG4P_E_);
0203             outtree->Branch("PrimaryG4P_PID", &PrimaryG4P_PID_);
0204             outtree->Branch("PrimaryG4P_trackID", &PrimaryG4P_trackID_);
0205             outtree->Branch("PrimaryG4P_isChargeHadron", &PrimaryG4P_isChargeHadron_);
0206             outtree->Branch("PHG4Hit_x0", &PHG4Hit_x0_);
0207             outtree->Branch("PHG4Hit_y0", &PHG4Hit_y0_);
0208             outtree->Branch("PHG4Hit_z0", &PHG4Hit_z0_);
0209             outtree->Branch("PHG4Hit_x1", &PHG4Hit_x1_);
0210             outtree->Branch("PHG4Hit_y1", &PHG4Hit_y1_);
0211             outtree->Branch("PHG4Hit_z1", &PHG4Hit_z1_);
0212             outtree->Branch("PHG4Hit_edep", &PHG4Hit_edep_);
0213             if (_get_allphg4_info)
0214             {
0215                 outtree->Branch("NAllG4P", &NAllG4P_);
0216                 outtree->Branch("G4P_Pt", &G4P_Pt_);
0217                 outtree->Branch("G4P_Eta", &G4P_Eta_);
0218                 outtree->Branch("G4P_Phi", &G4P_Phi_);
0219                 outtree->Branch("G4P_E", &G4P_E_);
0220                 outtree->Branch("G4P_PID", &G4P_PID_);
0221                 outtree->Branch("G4P_trackID", &G4P_trackID_);
0222             }
0223             // Truth cluster information & matching with reco clusters
0224             outtree->Branch("NTruthLayers", &NTruthLayers_);
0225             outtree->Branch("ClusTruthCKeys", &ClusTruthCKeys_);
0226             outtree->Branch("ClusNG4Particles", &ClusNG4Particles_);
0227             outtree->Branch("ClusNPrimaryG4Particles", &ClusNPrimaryG4Particles_);
0228             outtree->Branch("TruthClusPhiSize", &TruthClusPhiSize_);
0229             outtree->Branch("TruthClusZSize", &TruthClusZSize_);
0230             outtree->Branch("TruthClusNRecoClus", &TruthClusNRecoClus_);
0231             outtree->Branch("PrimaryTruthClusPhiSize", &PrimaryTruthClusPhiSize_);
0232             outtree->Branch("PrimaryTruthClusZSize", &PrimaryTruthClusZSize_);
0233             outtree->Branch("PrimaryTruthClusNRecoClus", &PrimaryTruthClusNRecoClus_);
0234             outtree->Branch("ClusMatchedG4P_MaxE_trackID", &ClusMatchedG4P_MaxE_trackID_);
0235             outtree->Branch("ClusMatchedG4P_MaxE_Pt", &ClusMatchedG4P_MaxE_Pt_);
0236             outtree->Branch("ClusMatchedG4P_MaxE_Eta", &ClusMatchedG4P_MaxE_Eta_);
0237             outtree->Branch("ClusMatchedG4P_MaxE_Phi", &ClusMatchedG4P_MaxE_Phi_);
0238             outtree->Branch("ClusMatchedG4P_MaxClusE_trackID", &ClusMatchedG4P_MaxClusE_trackID_);
0239             outtree->Branch("ClusMatchedG4P_MaxClusE_ancestorTrackID", &ClusMatchedG4P_MaxClusE_ancestorTrackID_);
0240             outtree->Branch("ClusMatchedG4P_MaxClusE_Pt", &ClusMatchedG4P_MaxClusE_Pt_);
0241             outtree->Branch("ClusMatchedG4P_MaxClusE_Eta", &ClusMatchedG4P_MaxClusE_Eta_);
0242             outtree->Branch("ClusMatchedG4P_MaxClusE_Phi", &ClusMatchedG4P_MaxClusE_Phi_);
0243         }
0244         // InttRawHit information
0245         if (_get_inttrawhit)
0246         {
0247             outtree->Branch("NInttRawHits", &NInttRawHits_);
0248             outtree->Branch("InttRawHit_bco", &InttRawHit_bco_);
0249             outtree->Branch("InttRawHit_packetid", &InttRawHit_packetid_);
0250             outtree->Branch("InttRawHit_word", &InttRawHit_word_);
0251             outtree->Branch("InttRawHit_fee", &InttRawHit_fee_);
0252             outtree->Branch("InttRawHit_channel_id", &InttRawHit_channel_id_);
0253             outtree->Branch("InttRawHit_chip_id", &InttRawHit_chip_id_);
0254             outtree->Branch("InttRawHit_adc", &InttRawHit_adc_);
0255             outtree->Branch("InttRawHit_FPHX_BCO", &InttRawHit_FPHX_BCO_);
0256             outtree->Branch("InttRawHit_full_FPHX", &InttRawHit_full_FPHX_);
0257             outtree->Branch("InttRawHit_full_ROC", &InttRawHit_full_ROC_);
0258             outtree->Branch("InttRawHit_amplitude", &InttRawHit_amplitude_);
0259         }
0260         // TrkrHit informatio
0261         if (_get_trkr_hit)
0262         {
0263             outtree->Branch("NTrkrhits", &NTrkrhits_);
0264             outtree->Branch("NTrkrhits_Layer1", &NTrkrhits_Layer1_);
0265             outtree->Branch("TrkrHitRow", &TrkrHitRow_);
0266             outtree->Branch("TrkrHitColumn", &TrkrHitColumn_);
0267             outtree->Branch("TrkrHitLadderZId", &TrkrHitLadderZId_);
0268             outtree->Branch("TrkrHitLadderPhiId", &TrkrHitLadderPhiId_);
0269             outtree->Branch("TrkrHitTimeBucketId", &TrkrHitTimeBucketId_);
0270             outtree->Branch("TrkrHitLayer", &TrkrHitLayer_);
0271             outtree->Branch("TrkrHitADC", &TrkrHitADC_);
0272             outtree->Branch("TrkrHitX", &TrkrHitX_);
0273             outtree->Branch("TrkrHitY", &TrkrHitY_);
0274             outtree->Branch("TrkrHitZ", &TrkrHitZ_);
0275             if (!IsData)
0276             {
0277                 // Truth PHG4Hit matching with TrkrHits
0278                 outtree->Branch("TrkrHit_truthHit_x0", &TrkrHit_truthHit_x0_);
0279                 outtree->Branch("TrkrHit_truthHit_y0", &TrkrHit_truthHit_y0_);
0280                 outtree->Branch("TrkrHit_truthHit_z0", &TrkrHit_truthHit_z0_);
0281                 outtree->Branch("TrkrHit_truthHit_x1", &TrkrHit_truthHit_x1_);
0282                 outtree->Branch("TrkrHit_truthHit_y1", &TrkrHit_truthHit_y1_);
0283                 outtree->Branch("TrkrHit_truthHit_z1", &TrkrHit_truthHit_z1_);
0284             }
0285         }
0286         // TrkrCluster information
0287         if (_get_reco_cluster)
0288         {
0289             outtree->Branch("NClus_Layer1", &NClus_Layer1_);
0290             outtree->Branch("NClus", &NClus_);
0291             outtree->Branch("ClusLayer", &ClusLayer_);
0292             outtree->Branch("ClusX", &ClusX_);
0293             outtree->Branch("ClusY", &ClusY_);
0294             outtree->Branch("ClusZ", &ClusZ_);
0295             outtree->Branch("ClusR", &ClusR_);
0296             outtree->Branch("ClusPhi", &ClusPhi_);
0297             outtree->Branch("ClusEta", &ClusEta_);
0298             outtree->Branch("ClusAdc", &ClusAdc_);
0299             outtree->Branch("ClusPhiSize", &ClusPhiSize_);
0300             outtree->Branch("ClusZSize", &ClusZSize_);
0301             outtree->Branch("ClusLadderZId", &ClusLadderZId_);
0302             outtree->Branch("ClusLadderPhiId", &ClusLadderPhiId_);
0303             outtree->Branch("ClusTrkrHitSetKey", &ClusTrkrHitSetKey_);
0304             outtree->Branch("ClusTimeBucketId", &ClusTimeBucketId_);
0305         }
0306         // GL1 trigger information
0307         if (_get_trigger_info)
0308         {
0309             outtree->Branch("GL1Packet_BCO", &GL1Packet_BCO_);
0310             outtree->Branch("triggervec", &triggervec_);
0311             outtree->Branch("firedTriggers", &firedTriggers_);
0312             outtree->Branch("firedTriggers_name", &firedTriggers_name_);
0313             outtree->Branch("firedTriggers_checkraw", &firedTriggers_checkraw_);
0314             outtree->Branch("firedTriggers_prescale", &firedTriggers_prescale_);
0315             outtree->Branch("firedTriggers_scalers", &firedTriggers_scalers_);
0316             outtree->Branch("firedTriggers_livescalers", &firedTriggers_livescalers_);
0317             outtree->Branch("firedTriggers_rawscalers", &firedTriggers_rawscalers_);
0318         }
0319     }
0320 
0321     return Fun4AllReturnCodes::EVENT_OK;
0322 }
0323 
0324 //____________________________________________________________________________..
0325 int dNdEtaINTT::InitRun(PHCompositeNode *topNode)
0326 {
0327     if (Verbosity() >= VERBOSITY_MORE)
0328         std::cout << "dNdEtaINTT::InitRun(PHCompositeNode *topNode) Initializing for Run XXX" << std::endl;
0329     return Fun4AllReturnCodes::EVENT_OK;
0330 }
0331 
0332 //____________________________________________________________________________..
0333 std::vector<int> dNdEtaINTT::GetAncestors(PHG4Particle *p)
0334 {
0335     // std::cout << "In GetG4Ancestor: " << std::endl << "start with original particle: ";
0336     // p->identify();
0337 
0338     // check if m_truth_info is available
0339     if (!m_truth_info)
0340     {
0341         std::cout << __func__ << " " << __LINE__ << ": PHG4TruthInfoContainer not found. Exit!" << std::endl;
0342         exit(1);
0343     }
0344 
0345     std::vector<int> ancestors_trackID;
0346     ancestors_trackID.clear();
0347 
0348     PHG4Particle *ancestor = m_truth_info->GetParticle(p->get_parent_id());
0349     // PHG4Particle *prevpart = p;
0350     while (ancestor != nullptr)
0351     {
0352         ancestors_trackID.push_back(ancestor->get_track_id());
0353         // ancestor->identify();
0354         // prevpart = ancestor;
0355         ancestor = m_truth_info->GetParticle(ancestor->get_parent_id());
0356     }
0357     return ancestors_trackID;
0358 }
0359 
0360 //____________________________________________________________________________..
0361 int dNdEtaINTT::process_event(PHCompositeNode *topNode)
0362 {
0363     if (Verbosity() >= VERBOSITY_MORE)
0364         std::cout << "dNdEtaINTT::process_event(PHCompositeNode *topNode) Processing Event " << eventNum << std::endl;
0365 
0366     PHNodeIterator nodeIter(topNode);
0367 
0368     if (!IsData)
0369     {
0370         if (!svtx_evalstack)
0371         {
0372             svtx_evalstack = new SvtxEvalStack(topNode);
0373             clustereval = svtx_evalstack->get_cluster_eval();
0374             hiteval = svtx_evalstack->get_hit_eval();
0375             truth_eval = svtx_evalstack->get_truth_eval();
0376         }
0377 
0378         svtx_evalstack->next_event(topNode);
0379     }
0380 
0381     if (_get_centrality)
0382     {
0383         eventheader = findNode::getClass<EventHeader>(topNode, "EventHeader");
0384         if (!eventheader)
0385         {
0386             std::cout << "Error, can't find EventHeader" << std::endl;
0387             exit(1);
0388         }
0389 
0390         GetCentralityInfo(topNode);
0391     }
0392 
0393     if (_get_intt_data)
0394     {
0395         if (!IsData)
0396         {
0397             std::cout << "Running on simulation" << std::endl;
0398 
0399             if (_get_phg4_info)
0400                 GetPHG4Info(topNode);
0401 
0402             if (_get_allphg4_info)
0403                 GetAllPHG4Info(topNode);
0404 
0405             if (_get_hepmc_info)
0406             {
0407                 // std::cout << "[INFO & WARNING] Currently, GetHEPMCInfo method is disabled." << std::endl;
0408                 GetHEPMCInfo(topNode);
0409             }
0410 
0411             if (_get_truth_cluster)
0412                 GetTruthClusterInfo(topNode);
0413         }
0414 
0415         if (_get_inttrawhit)
0416             GetInttRawHitInfo(topNode);
0417 
0418         if (_get_trkr_hit)
0419             GetTrkrHitInfo(topNode);
0420 
0421         if (_get_reco_cluster)
0422             GetRecoClusterInfo(topNode);
0423     }
0424 
0425     if (IsData)
0426     {
0427         if (_get_intt_data)
0428         {
0429             intteventinfo = findNode::getClass<InttEventInfo>(topNode, "INTTEVENTHEADER");
0430             if (!intteventinfo)
0431             {
0432                 std::cout << "The INTT event header is missing, you will have no BCO information fro syncing" << std::endl;
0433             }
0434 
0435             intt_bco = intteventinfo->get_bco_full();
0436         }
0437 
0438         if (_get_trigger_info)
0439         {
0440             gl1packet = findNode::getClass<Gl1Packet>(topNode, "GL1RAWHIT");
0441             if (!gl1packet)
0442             {
0443                 std::cout << __FILE__ << "::" << __func__ << " - Gl1Packet missing, doing nothing." << std::endl;
0444                 exit(1);
0445             }
0446             GetTriggerInfo(topNode);
0447         }
0448     }
0449 
0450     // event_ = InputFileListIndex * NEvtPerFile + eventNum;
0451     event_ = eventNum;
0452     // evt_sequence_ = (_get_centrality) ? eventheader->get_EvtSequence() : eventNum;
0453     outtree->Fill();
0454     eventNum++;
0455 
0456     ResetVectors();
0457 
0458     return Fun4AllReturnCodes::EVENT_OK;
0459 }
0460 
0461 //____________________________________________________________________________..
0462 int dNdEtaINTT::ResetEvent(PHCompositeNode *topNode)
0463 {
0464     if (Verbosity() >= VERBOSITY_MORE)
0465         std::cout << "dNdEtaINTT::ResetEvent(PHCompositeNode *topNode) Resetting internal structures, prepare for next event" << std::endl;
0466     return Fun4AllReturnCodes::EVENT_OK;
0467 }
0468 
0469 //____________________________________________________________________________..
0470 int dNdEtaINTT::EndRun(const int runnumber)
0471 {
0472     if (Verbosity() >= VERBOSITY_MORE)
0473         std::cout << "dNdEtaINTT::EndRun(const int runnumber) Ending Run for Run " << runnumber << std::endl;
0474     return Fun4AllReturnCodes::EVENT_OK;
0475 }
0476 
0477 //____________________________________________________________________________..
0478 int dNdEtaINTT::End(PHCompositeNode *topNode)
0479 {
0480     if (Verbosity() >= VERBOSITY_MORE)
0481         std::cout << "dNdEtaINTT::End(PHCompositeNode *topNode) This is the End - Output to " << _outputFile << std::endl;
0482 
0483     PHTFileServer::get().cd(_outputFile);
0484     outtree->Write("", TObject::kOverwrite);
0485 
0486     delete svtx_evalstack;
0487 
0488     return Fun4AllReturnCodes::EVENT_OK;
0489 }
0490 
0491 //____________________________________________________________________________..
0492 int dNdEtaINTT::Reset(PHCompositeNode *topNode)
0493 {
0494     if (Verbosity() >= VERBOSITY_MORE)
0495         std::cout << "dNdEtaINTT::Reset(PHCompositeNode *topNode) being Reset" << std::endl;
0496     return Fun4AllReturnCodes::EVENT_OK;
0497 }
0498 
0499 //____________________________________________________________________________..
0500 void dNdEtaINTT::Print(const std::string &what) const { std::cout << "dNdEtaINTT::Print(const std::string &what) const Printing info for " << what << std::endl; }
0501 //____________________________________________________________________________..
0502 void dNdEtaINTT::GetTriggerInfo(PHCompositeNode *topNode)
0503 {
0504     std::cout << "Get GL1 trigger info." << std::endl;
0505 
0506     GL1Packet_BCO_ = (gl1packet->lValue(0, "BCO") & 0xFFFFFFFFFFU); // Only the lower 40 bits are retained
0507 
0508     firedTriggers_.clear();
0509     firedTriggers_name_.clear();
0510     firedTriggers_checkraw_.clear();
0511     firedTriggers_prescale_.clear();
0512     firedTriggers_scalers_.clear();
0513     firedTriggers_livescalers_.clear();
0514     firedTriggers_rawscalers_.clear();
0515 
0516     // std::cout << __FILE__ << "::" << __func__ << "::" << __LINE__ << std::endl;
0517 
0518     // Set up Trigger Analyzer
0519     triggeranalyzer = new TriggerAnalyzer();
0520     triggeranalyzer->decodeTriggers(topNode);
0521 
0522     // std::cout << __FILE__ << "::" << __func__ << "::" << __LINE__ << std::endl;
0523 
0524     // triggervec_ = gl1packet->getTriggerVector(); // just to get the original triggervec
0525     triggervec_ = gl1packet->getScaledVector();
0526 
0527     // std::cout << __FILE__ << "::" << __func__ << "::" << __LINE__ << std::endl;
0528 
0529     for (int i = 0; i < 64; i++)
0530     {
0531         bool trig_decision = ((triggervec_ & 0x1U) == 0x1U);
0532         if (trig_decision)
0533         {
0534             firedTriggers_.push_back(i);
0535             firedTriggers_name_.push_back(triggeranalyzer->getTriggerName(i));
0536             firedTriggers_checkraw_.push_back(triggeranalyzer->checkRawTrigger(i));
0537             firedTriggers_prescale_.push_back(triggeranalyzer->getTriggerPrescale(i));
0538             firedTriggers_scalers_.push_back(triggeranalyzer->getTriggerScalers(i));
0539             firedTriggers_livescalers_.push_back(triggeranalyzer->getTriggerLiveScalers(i));
0540             firedTriggers_rawscalers_.push_back(triggeranalyzer->getTriggerRawScalers(i));
0541         }
0542         triggervec_ = (triggervec_ >> 1U) & 0xffffffffU;
0543     }
0544 
0545     // std::cout << __FILE__ << "::" << __func__ << "::" << __LINE__ << std::endl;
0546 
0547     triggervec_ = gl1packet->getScaledVector(); // just to get the original triggervec
0548 }
0549 //____________________________________________________________________________..
0550 //! error when compiling in ALMA9 with c++20 -> change -IOFFLINEMAIN/include to -isystemOFFLINEMAIN/include in Makefile as a workaround
0551 void dNdEtaINTT::GetHEPMCInfo(PHCompositeNode *topNode)
0552 {
0553     std::cout << "Get HEPMC info." << std::endl;
0554 
0555     // HEPMC info
0556     m_geneventmap = findNode::getClass<PHHepMCGenEventMap>(topNode, "PHHepMCGenEventMap");
0557     PHHepMCGenHelper *phhepmcgenhlpr = new PHHepMCGenHelper();
0558     phhepmcgenhlpr->PHHepMCGenHelper_Verbosity(1);
0559 
0560     if (m_geneventmap)
0561     {
0562         PHHepMCGenEventMap::ConstIter iter_beg = m_geneventmap->begin();
0563         PHHepMCGenEventMap::ConstIter iter_end = m_geneventmap->end();
0564         for (PHHepMCGenEventMap::ConstIter iter = iter_beg; iter != iter_end; ++iter)
0565         {
0566             PHHepMCGenEvent *genevt = iter->second;
0567             // HepMC::ThreeVector initv_boostbeta = genevt->get_boost_beta_vector(); HepMC::ThreeVector
0568             // initv_rotation = genevt->get_rotation_vector();
0569             // std::cout << "Initial boost beta vector: " << initv_boostbeta.x() << " " << initv_boostbeta.y() << " " << initv_boostbeta.z() << std::endl;
0570             // std::cout << "Initial rotation vector: " << initv_rotation.x() << " " << initv_rotation.y() << " " << initv_rotation.z() << std::endl;
0571             // genevt->identify();
0572 
0573             if (genevt->get_embedding_id() != 0)
0574                 continue;
0575 
0576             HepMC::GenEvent *evt = genevt->getEvent();
0577 
0578             // Get the Lorentz rotation and boost
0579             const CLHEP::HepLorentzRotation lortentz_rotation(genevt->get_LorentzRotation_EvtGen2Lab());
0580             const double mom_factor = HepMC::Units::conversion_factor(evt->momentum_unit(), HepMC::Units::GEV);
0581 
0582             if (evt)
0583             {
0584                 // genevt->PrintEvent();
0585                 std::cout << "Signal process id " << evt->signal_process_id() << std::endl;
0586                 std::cout << "Embedding id " << genevt->get_embedding_id() << std::endl;
0587                 std::cout << "is simulated " << genevt->is_simulated() << std::endl;
0588                 std::cout << "Collision vertex x (PHHepMCGenEvent): " << genevt->get_collision_vertex().x() << std::endl;
0589                 std::cout << "Collision vertex y (PHHepMCGenEvent): " << genevt->get_collision_vertex().y() << std::endl;
0590                 std::cout << "Collision vertex z (PHHepMCGenEvent): " << genevt->get_collision_vertex().z() << std::endl;
0591                 std::cout << "Collision vertex t (PHHepMCGenEvent): " << genevt->get_collision_vertex().t() << std::endl;
0592 
0593                 signal_process_id_ = evt->signal_process_id();
0594 
0595                 for (HepMC::GenEvent::particle_const_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p)
0596                 {
0597                     HepMC::GenParticle *particle = *p;
0598                     HepMC::GenVertex *vtx_decay = particle->end_vertex();
0599                     if (vtx_decay == nullptr && particle->status() == 1)
0600                     {
0601                         // Lorentz rotation and boost for the beam crossing and shifted vertex
0602                         CLHEP::HepLorentzVector lv_momentum(particle->momentum().px(), particle->momentum().py(), particle->momentum().pz(), particle->momentum().e());
0603                         lv_momentum = lortentz_rotation(lv_momentum);
0604 
0605                         TLorentzVector hepmcp;
0606                         hepmcp.SetPxPyPzE(lv_momentum.px() * mom_factor, lv_momentum.py() * mom_factor, lv_momentum.pz() * mom_factor, lv_momentum.e() * mom_factor);
0607 
0608                         HepMCFSPrtl_Pt_.push_back(hepmcp.Pt());
0609                         HepMCFSPrtl_Eta_.push_back(hepmcp.Eta());
0610                         HepMCFSPrtl_Phi_.push_back(hepmcp.Phi());
0611                         HepMCFSPrtl_E_.push_back(hepmcp.E());
0612                         HepMCFSPrtl_PID_.push_back(particle->pdg_id());
0613 
0614                         HepMC::GenVertex *vtx_prod = particle->production_vertex();
0615                         if (vtx_prod)
0616                         {
0617                             HepMCFSPrtl_prodx_.push_back(vtx_prod->position().x());
0618                             HepMCFSPrtl_prody_.push_back(vtx_prod->position().y());
0619                             HepMCFSPrtl_prodz_.push_back(vtx_prod->position().z());
0620                         }
0621                     }
0622                 }
0623 
0624                 NHepMCFSPart_ = HepMCFSPrtl_PID_.size();
0625             }
0626             else
0627             {
0628                 std::cout << "No HepMC event" << std::endl;
0629                 continue;
0630             }
0631         }
0632     }
0633     else
0634     {
0635         std::cout << PHWHERE << "Error, can't find PHHepMCGenEventMap" << std::endl;
0636         exit(1);
0637     }
0638 }
0639 
0640 //____________________________________________________________________________..
0641 void dNdEtaINTT::GetCentralityInfo(PHCompositeNode *topNode)
0642 {
0643     if (Verbosity() >= VERBOSITY_MORE)
0644         std::cout << "Get centrality info." << std::endl;
0645 
0646     m_CentInfo = findNode::getClass<CentralityInfo>(topNode, "CentralityInfo");
0647     if (!m_CentInfo)
0648     {
0649         std::cout << PHWHERE << "Error, can't find CentralityInfov1. No centrality info is filled" << std::endl;
0650         // exit(1);
0651     }
0652 
0653     _minimumbiasinfo = findNode::getClass<MinimumBiasInfo>(topNode, "MinimumBiasInfo");
0654     if (!_minimumbiasinfo)
0655     {
0656         std::cout << "Error, can't find MinimumBiasInfo. No minimum bias info is filled" << std::endl;
0657         // exit(1);
0658     }
0659 
0660     m_mbdout = findNode::getClass<MbdOut>(topNode, "MbdOut");
0661     if (!m_mbdout)
0662     {
0663         std::cout << "Error, can't find MbdOut" << std::endl;
0664         exit(1);
0665     }
0666 
0667     if (_get_pmt_info)
0668     {
0669         m_mbdpmtcontainer = findNode::getClass<MbdPmtContainer>(topNode, "MbdPmtContainer");
0670         if (!m_mbdpmtcontainer)
0671         {
0672             std::cout << "Error, can't find MbdPmtContainer" << std::endl;
0673             exit(1);
0674         }
0675     }
0676 
0677     m_mbdvtxmap = findNode::getClass<MbdVertexMapv1>(topNode, "MbdVertexMap");
0678     if (!m_mbdvtxmap)
0679     {
0680         std::cout << "Error, can't find MbdVertexMap" << std::endl;
0681         exit(1);
0682     }
0683 
0684     mbd_z_vtx = std::numeric_limits<float>::quiet_NaN();
0685     std::cout << "MbdVertexMap size: " << m_mbdvtxmap->size() << std::endl;
0686     for (MbdVertexMap::ConstIter biter = m_mbdvtxmap->begin(); biter != m_mbdvtxmap->end(); ++biter)
0687     {
0688         m_mbdvtx = biter->second;
0689         mbd_z_vtx = m_mbdvtx->get_z();
0690     }
0691 
0692     if (!IsData)
0693     {
0694         centrality_bimp_ = (m_CentInfo) ? m_CentInfo->get_centile(CentralityInfo::PROP::bimp) : std::numeric_limits<float>::quiet_NaN();
0695         centrality_impactparam_ = (m_CentInfo) ? m_CentInfo->get_quantity(CentralityInfo::PROP::bimp) : std::numeric_limits<float>::quiet_NaN();
0696         // Glauber parameter information
0697         ncoll_ = eventheader->get_ncoll();
0698         npart_ = eventheader->get_npart();
0699         std::cout << "Centrality: (bimp,impactparam) = (" << centrality_bimp_ << ", " << centrality_impactparam_ << ")" << std::endl;
0700         std::cout << "Glauber parameter information: (ncoll,npart) = (" << ncoll_ << ", " << npart_ << ")" << std::endl;
0701     }
0702 
0703     clk = (IsData) ? m_mbdout->get_clock() : 0;
0704     femclk = (IsData) ? m_mbdout->get_femclock() : 0;
0705     mbd_south_npmt = m_mbdout->get_npmt(0);
0706     mbd_north_npmt = m_mbdout->get_npmt(1);
0707     mbd_south_charge_sum = m_mbdout->get_q(0);
0708     mbd_north_charge_sum = m_mbdout->get_q(1);
0709     mbd_charge_sum = mbd_south_charge_sum + mbd_north_charge_sum;
0710     mbd_charge_asymm = mbd_charge_sum == 0 ? std::numeric_limits<float>::quiet_NaN() : (float)(mbd_south_charge_sum - mbd_north_charge_sum) / mbd_charge_sum;
0711     if (m_CentInfo)
0712     {
0713         if (m_CentInfo->has_centrality_bin(CentralityInfo::PROP::mbd_NS))
0714         {
0715             centrality_mbd_ = m_CentInfo->get_centrality_bin(CentralityInfo::PROP::mbd_NS);
0716         }
0717         else
0718         {
0719             std::cout << "[WARNING/ERROR] No centrality information found in CentralityInfo. Setting centrality_mbd_ to NaN. Please check!" << std::endl;
0720             m_CentInfo->identify();
0721             centrality_mbd_ = std::numeric_limits<float>::quiet_NaN();
0722         }
0723     }
0724     else
0725     {
0726         centrality_mbd_ = std::numeric_limits<float>::quiet_NaN();
0727     }
0728 
0729     is_min_bias = (_minimumbiasinfo) ? _minimumbiasinfo->isAuAuMinimumBias() : false;
0730 
0731     std::cout << "[INFO] Is minimum bias: " << is_min_bias << "; Centrality: (centrality_mbd, centrality_mbdquantity) = (" << centrality_mbd_ << ", " << centrality_mbdquantity_ << ")" << std::endl;
0732 
0733     // minimum bias criteria without zdc cut (note: the zdc cut has a 99-100% efficiency for the Level-1 trigger events and 100% for central Au+Au event)
0734     bool mbd_ntube = (mbd_south_npmt >= 2 && mbd_north_npmt >= 2) ? true : false;
0735     bool mbd_sn_q_imbalence = (mbd_north_charge_sum > 10 || mbd_south_charge_sum < 150) ? true : false;
0736     bool mbd_zvtx = (fabs(mbd_z_vtx) < 60) ? true : false;
0737     is_min_bias_wozdc = (IsData) ? (mbd_ntube && mbd_sn_q_imbalence && mbd_zvtx) : (npart_ > 0);
0738 
0739     int hits_n = 0;
0740     int hits_s = 0;
0741 
0742     if (_get_pmt_info)
0743     {
0744         for (int i = 0; i < 128; i++)
0745         {
0746             m_pmt_q[i] = m_mbdpmtcontainer->get_pmt(i)->get_q();
0747             if (i < 64 && m_pmt_q[i] > cthresh)
0748                 hits_s++;
0749             else if (i >= 64 && m_pmt_q[i] > cthresh)
0750                 hits_n++;
0751         }
0752     }
0753     mbd_nhitsoverths_south = hits_s;
0754     mbd_nhitsoverths_north = hits_n;
0755 }
0756 //____________________________________________________________________________..
0757 void dNdEtaINTT::GetInttRawHitInfo(PHCompositeNode *topNode)
0758 {
0759     std::cout << "Get InttRawHit info." << std::endl;
0760 
0761     inttrawhitcontainer = findNode::getClass<InttRawHitContainer>(topNode, "INTTRAWHIT");
0762     if (!inttrawhitcontainer)
0763     {
0764         std::cout << PHWHERE << "Error, can't find INTTRAWHIT" << std::endl;
0765         exit(1);
0766     }
0767 
0768     NInttRawHits_ = inttrawhitcontainer->get_nhits();
0769 
0770     for (unsigned int i = 0; i < inttrawhitcontainer->get_nhits(); i++)
0771     {
0772         InttRawHit *intthit = inttrawhitcontainer->get_hit(i);
0773 
0774         InttRawHit_bco_.push_back(intthit->get_bco());
0775         InttRawHit_packetid_.push_back(intthit->get_packetid());
0776         InttRawHit_word_.push_back(intthit->get_word());
0777         InttRawHit_fee_.push_back(intthit->get_fee());
0778         InttRawHit_channel_id_.push_back(intthit->get_channel_id());
0779         InttRawHit_chip_id_.push_back(intthit->get_chip_id());
0780         InttRawHit_adc_.push_back(intthit->get_adc());
0781         InttRawHit_FPHX_BCO_.push_back(intthit->get_FPHX_BCO());
0782         InttRawHit_full_FPHX_.push_back(intthit->get_full_FPHX());
0783         InttRawHit_full_ROC_.push_back(intthit->get_full_ROC());
0784         InttRawHit_amplitude_.push_back(intthit->get_amplitude());
0785     }
0786 }
0787 //____________________________________________________________________________..
0788 void dNdEtaINTT::GetTrkrHitInfo(PHCompositeNode *topNode)
0789 {
0790     if (Verbosity() >= VERBOSITY_MORE)
0791         std::cout << "Get TrkrHit info." << std::endl;
0792 
0793     hitsets = findNode::getClass<TrkrHitSetContainer>(topNode, "TRKR_HITSET");
0794     if (!hitsets)
0795     {
0796         std::cout << PHWHERE << "Error, can't find cluster-hit map" << std::endl;
0797         exit(1);
0798     }
0799 
0800     _intt_geom_container = findNode::getClass<PHG4CylinderGeomContainer>(topNode, "CYLINDERGEOM_INTT");
0801     if (!_intt_geom_container)
0802     {
0803         std::cout << PHWHERE << "Error, can't find INTT geometry container" << std::endl;
0804         exit(1);
0805     }
0806 
0807     _tgeometry = findNode::getClass<ActsGeometry>(topNode, "ActsGeometry");
0808     if (!_tgeometry)
0809     {
0810         std::cout << PHWHERE << "Error, can't find Acts geometry" << std::endl;
0811         exit(1);
0812     }
0813 
0814     if (!IsData)
0815     {
0816         _hit_truth_map = findNode::getClass<TrkrHitTruthAssoc>(topNode, "TRKR_HITTRUTHASSOC");
0817         if (!_hit_truth_map)
0818         {
0819             std::cout << PHWHERE << "Error, can't find TRKR_HITTRUTHASSOC" << std::endl;
0820             exit(1);
0821         }
0822 
0823         g4hit = findNode::getClass<PHG4HitContainer>(topNode, "G4HIT_INTT");
0824         if (!g4hit)
0825         {
0826             std::cout << PHWHERE << "Error, can't find G4HIT_INTT" << std::endl;
0827             exit(1);
0828         }
0829     }
0830 
0831     std::set<PHG4Hit *> truth_hits; // unique set of matched PHG4Hits
0832     // make sure the set is empty
0833     truth_hits.clear();
0834 
0835     TrkrHitSetContainer::ConstRange hitset_range = hitsets->getHitSets(TrkrDefs::TrkrId::inttId);
0836     for (TrkrHitSetContainer::ConstIterator hitset_iter = hitset_range.first; hitset_iter != hitset_range.second; ++hitset_iter)
0837     {
0838         TrkrHitSet::ConstRange hit_range = hitset_iter->second->getHits();
0839         TrkrHitSet *hitset = hitset_iter->second;
0840         auto hitsetkey = hitset->getHitSetKey();
0841 
0842         auto surface = _tgeometry->maps().getSiliconSurface(hitsetkey);
0843         auto layergeom = dynamic_cast<CylinderGeomIntt *>(_intt_geom_container->GetLayerGeom(TrkrDefs::getLayer(hitsetkey)));
0844         TVector2 LocalUse;
0845 
0846         for (TrkrHitSet::ConstIterator hit_iter = hit_range.first; hit_iter != hit_range.second; ++hit_iter)
0847         {
0848             TrkrDefs::hitkey hitKey = hit_iter->first;
0849 
0850             // now get the positions from the geometry
0851             auto col = InttDefs::getCol(hitKey);
0852             auto row = InttDefs::getRow(hitKey);
0853             auto ladder_z_index = InttDefs::getLadderZId(hitsetkey);
0854             double local_hit_location[3] = {0., 0., 0.};
0855             layergeom->find_strip_center_localcoords(ladder_z_index, row, col, local_hit_location);
0856             LocalUse.SetX(local_hit_location[0]); // local x
0857             LocalUse.SetY(local_hit_location[2]); // local y
0858             TVector3 posworld = CylinderGeomInttHelper::get_world_from_local_coords(surface, _tgeometry, LocalUse);
0859 
0860             // std::cout << "Hit position: " << posworld.x() << " " << posworld.y() << " " << posworld.z() << std::endl;
0861 
0862             TrkrHitX_.push_back(posworld.x());
0863             TrkrHitY_.push_back(posworld.y());
0864             TrkrHitZ_.push_back(posworld.z());
0865             TrkrHitRow_.push_back(InttDefs::getRow(hitKey));
0866             TrkrHitColumn_.push_back(InttDefs::getCol(hitKey));
0867             TrkrHitLadderZId_.push_back(InttDefs::getLadderZId(hitsetkey));
0868             TrkrHitLadderPhiId_.push_back(InttDefs::getLadderPhiId(hitsetkey));
0869             TrkrHitTimeBucketId_.push_back(InttDefs::getTimeBucketId(hitsetkey));
0870             TrkrHitLayer_.push_back(TrkrDefs::getLayer(hitsetkey));
0871             TrkrHitADC_.push_back(hit_iter->second->getAdc());
0872 
0873             // https://github.com/sPHENIX-Collaboration/coresoftware/blob/master/simulation/g4simulation/g4eval/SvtxHitEval.cc
0874             // std::set<PHG4Hit*> truth_hits;
0875             if (!IsData)
0876             {
0877                 TrkrHit *hit = hit_iter->second;
0878                 if (hit)
0879                 {
0880                     std::multimap<TrkrDefs::hitsetkey, std::pair<TrkrDefs::hitkey, PHG4HitDefs::keytype>> temp_map;
0881                     _hit_truth_map->getG4Hits(hitsetkey, hitKey, temp_map); // returns pairs (hitsetkey, std::pair(hitkey, g4hitkey)) for this hitkey only
0882                     for (auto &htiter : temp_map)
0883                     {
0884                         // extract the g4 hit key here and add the g4hit to the set
0885                         PHG4HitDefs::keytype g4hitkey = htiter.second.second;
0886                         // std::cout << "           hitkey " << hitkey <<  " g4hitkey " << g4hitkey << std::endl;
0887                         PHG4Hit *trkrhit_truthhit = nullptr;
0888                         trkrhit_truthhit = g4hit->findHit(g4hitkey);
0889                         if (trkrhit_truthhit)
0890                         {
0891                             truth_hits.insert(trkrhit_truthhit);
0892                         }
0893                     }
0894                 }
0895             }
0896         }
0897     }
0898 
0899     NTrkrhits_ = TrkrHitRow_.size();
0900     NTrkrhits_Layer1_ = std::count_if(TrkrHitLayer_.begin(), TrkrHitLayer_.end(), [](int i) { return i == 3 || i == 4; });
0901 
0902     if (!IsData)
0903     {
0904         for (auto &hit : truth_hits)
0905         {
0906             TrkrHit_truthHit_x0_.push_back(hit->get_x(0));
0907             TrkrHit_truthHit_y0_.push_back(hit->get_y(0));
0908             TrkrHit_truthHit_z0_.push_back(hit->get_z(0));
0909             TrkrHit_truthHit_x1_.push_back(hit->get_x(1));
0910             TrkrHit_truthHit_y1_.push_back(hit->get_y(1));
0911             TrkrHit_truthHit_z1_.push_back(hit->get_z(1));
0912         }
0913     }
0914 }
0915 //____________________________________________________________________________..
0916 void dNdEtaINTT::GetRecoClusterInfo(PHCompositeNode *topNode)
0917 {
0918     if (Verbosity() >= VERBOSITY_MORE)
0919         std::cout << "Get reconstructed cluster info." << std::endl;
0920 
0921     dst_clustermap = findNode::getClass<TrkrClusterContainerv4>(topNode, "TRKR_CLUSTER");
0922     if (!dst_clustermap)
0923     {
0924         std::cout << PHWHERE << "Error, can't find trkr clusters" << std::endl;
0925         exit(1);
0926     }
0927 
0928     _tgeometry = findNode::getClass<ActsGeometry>(topNode, "ActsGeometry");
0929     if (!_tgeometry)
0930     {
0931         std::cout << PHWHERE << "Error, can't find Acts geometry" << std::endl;
0932         exit(1);
0933     }
0934 
0935     std::vector<int> _NClus = {0, 0};
0936     // std::map<int, unsigned int> AncG4P_cluster_count;
0937     // AncG4P_cluster_count.clear();
0938 
0939     // Reference:
0940     // https://github.com/sPHENIX-Collaboration/coresoftware/blob/master/simulation/g4simulation/g4eval/SvtxEvaluator.cc#L1731-L1984
0941     // for (const auto &hitsetkey : dst_clustermap->getHitSetKeys())
0942     for (const auto &hitsetkey : dst_clustermap->getHitSetKeys(TrkrDefs::inttId))
0943     {
0944         auto range = dst_clustermap->getClusters(hitsetkey);
0945         for (auto iter = range.first; iter != range.second; ++iter)
0946         {
0947             // std::cout << "----------" << std::endl;
0948             TrkrDefs::cluskey ckey = iter->first;
0949             TrkrCluster *cluster = dst_clustermap->findCluster(ckey);
0950             unsigned int trkrId = TrkrDefs::getTrkrId(ckey);
0951             if (trkrId != TrkrDefs::inttId)
0952                 continue; // we want only INTT clusters
0953             int layer = (TrkrDefs::getLayer(ckey) == 3 || TrkrDefs::getLayer(ckey) == 4) ? 0 : 1;
0954             _NClus[layer]++;
0955             if (cluster == nullptr)
0956             {
0957                 std::cout << "cluster is nullptr, ckey=" << ckey << std::endl;
0958             }
0959             auto globalpos = _tgeometry->getGlobalPosition(ckey, cluster);
0960             ClusLayer_.push_back(TrkrDefs::getLayer(ckey));
0961             ClusX_.push_back(globalpos(0));
0962             ClusY_.push_back(globalpos(1));
0963             ClusZ_.push_back(globalpos(2));
0964             ClusAdc_.push_back(cluster->getAdc());
0965             TVector3 pos(globalpos(0), globalpos(1), globalpos(2));
0966             ClusR_.push_back(pos.Perp());
0967             ClusPhi_.push_back(pos.Phi());
0968             ClusEta_.push_back(pos.Eta());
0969             // ClusPhiSize is a signed char (-127-127), we should convert it to an unsigned char (0-255)
0970             // ClusPhiSize is a signed char (-127-127), we should convert it to
0971             // an unsigned char (0-255)
0972             int phisize = cluster->getPhiSize();
0973             if (phisize <= 0)
0974                 phisize += 256;
0975             ClusPhiSize_.push_back(phisize);
0976             ClusZSize_.push_back(cluster->getZSize());
0977             ClusLadderZId_.push_back(InttDefs::getLadderZId(ckey));
0978             ClusLadderPhiId_.push_back(InttDefs::getLadderPhiId(ckey));
0979             ClusTrkrHitSetKey_.push_back(hitsetkey);
0980             ClusTimeBucketId_.push_back(InttDefs::getTimeBucketId(ckey));
0981             if (!IsData)
0982             {
0983                 // truth cluster association
0984                 std::vector<int32_t> truth_ckeys;
0985                 std::set<PHG4Particle *> truth_particles = clustereval->all_truth_particles(ckey);
0986                 int Nprimary = 0;
0987                 for (auto &p : truth_particles)
0988                 {
0989                     // must be primary truth particle
0990                     if (p->get_parent_id() == 0)
0991                     {
0992                         Nprimary++;
0993                         std::map<TrkrDefs::cluskey, std::shared_ptr<TrkrCluster>> truth_clusters = truth_eval->all_truth_clusters(p);
0994                         for (auto &c : truth_clusters)
0995                         {
0996                             // only take at most 1 truth cluster per truth particle
0997                             bool found = false;
0998                             if (TrkrDefs::getLayer(c.first) == TrkrDefs::getLayer(ckey))
0999                             {
1000                                 if (!found)
1001                                 {
1002                                     found = true;
1003                                     truth_ckeys.push_back(c.first);
1004                                 }
1005                             }
1006                         }
1007                     }
1008                 }
1009                 ClusTruthCKeys_.push_back(truth_ckeys.size());
1010                 ClusNG4Particles_.push_back(truth_particles.size());
1011                 ClusNPrimaryG4Particles_.push_back(Nprimary);
1012 
1013                 // max_truth_cluster_by_energy
1014                 PHG4Particle *ptcl_maxE = clustereval->max_truth_particle_by_energy(ckey);
1015                 if (ptcl_maxE)
1016                 {
1017                     ClusMatchedG4P_MaxE_trackID_.push_back(ptcl_maxE->get_track_id());
1018                     TLorentzVector p;
1019                     p.SetPxPyPzE(ptcl_maxE->get_px(), ptcl_maxE->get_py(), ptcl_maxE->get_pz(), ptcl_maxE->get_e());
1020                     ClusMatchedG4P_MaxE_Pt_.push_back(p.Pt());
1021                     ClusMatchedG4P_MaxE_Eta_.push_back(p.Eta());
1022                     ClusMatchedG4P_MaxE_Phi_.push_back(p.Phi());
1023                 }
1024                 else
1025                 {
1026                     ClusMatchedG4P_MaxE_trackID_.push_back(std::numeric_limits<int>::max());
1027                     ClusMatchedG4P_MaxE_Pt_.push_back(-1);
1028                     ClusMatchedG4P_MaxE_Eta_.push_back(-1);
1029                     ClusMatchedG4P_MaxE_Phi_.push_back(-1);
1030                 }
1031                 // max_truth_particle_by_cluster_energy
1032                 PHG4Particle *ptcl_maxClusE = clustereval->max_truth_particle_by_cluster_energy(ckey);
1033                 if (ptcl_maxClusE)
1034                 {
1035                     ClusMatchedG4P_MaxClusE_trackID_.push_back(ptcl_maxClusE->get_track_id());
1036                     TLorentzVector p;
1037                     p.SetPxPyPzE(ptcl_maxClusE->get_px(), ptcl_maxClusE->get_py(), ptcl_maxClusE->get_pz(), ptcl_maxClusE->get_e());
1038                     ClusMatchedG4P_MaxClusE_Pt_.push_back(p.Pt());
1039                     ClusMatchedG4P_MaxClusE_Eta_.push_back(p.Eta());
1040                     ClusMatchedG4P_MaxClusE_Phi_.push_back(p.Phi());
1041                     // ancestors of the max cluster energy truth particle
1042                     std::vector<int> ptcl_maxClusE_ancestorsTrackID = GetAncestors(ptcl_maxClusE);
1043                     // if no ancestors, then the trackID of the max cluster energy truth particle is the ancestor
1044                     // if has ancestors, then the trackID of the last element of ptcl_maxClusE_ancestorsTrackID is the ancestor trackID
1045                     ClusMatchedG4P_MaxClusE_ancestorTrackID_.push_back((ptcl_maxClusE_ancestorsTrackID.size() == 0) ? ptcl_maxClusE->get_track_id() : ptcl_maxClusE_ancestorsTrackID.back());
1046                     // std::cout << "Macthed ptcl_maxClusE trackID: " << ptcl_maxClusE->get_track_id() << " - Ancestors of the max cluster energy truth particle: ";
1047                     // for (auto &i : ptcl_maxClusE_ancestorsTrackID)
1048                     //     std::cout << i << " ";
1049                     // std::cout << std::endl;
1050                 }
1051                 else
1052                 {
1053                     ClusMatchedG4P_MaxClusE_trackID_.push_back(std::numeric_limits<int>::max());
1054                     ClusMatchedG4P_MaxClusE_Pt_.push_back(-1);
1055                     ClusMatchedG4P_MaxClusE_Eta_.push_back(-1);
1056                     ClusMatchedG4P_MaxClusE_Phi_.push_back(-1);
1057                 }
1058             }
1059         }
1060     }
1061 
1062     NClus_ = _NClus[0] + _NClus[1];
1063     NClus_Layer1_ = _NClus[0];
1064     if (Verbosity() >= VERBOSITY_MORE)
1065         std::cout << "Number of clusters (total,0,1)=(" << NClus_ << "," << _NClus[0] << "," << _NClus[1] << ")" << std::endl;
1066 }
1067 //____________________________________________________________________________..
1068 void dNdEtaINTT::GetTruthClusterInfo(PHCompositeNode *topNode)
1069 {
1070     std::cout << "Get truth cluster info." << std::endl;
1071 
1072     _intt_geom_container = findNode::getClass<PHG4CylinderGeomContainer>(topNode, "CYLINDERGEOM_INTT");
1073     if (!_intt_geom_container)
1074     {
1075         std::cout << PHWHERE << "Error, can't find INTT geometry container" << std::endl;
1076         exit(1);
1077     }
1078 
1079     _tgeometry = findNode::getClass<ActsGeometry>(topNode, "ActsGeometry");
1080     if (!_tgeometry)
1081     {
1082         std::cout << PHWHERE << "Error, can't find Acts geometry" << std::endl;
1083         exit(1);
1084     }
1085 
1086     auto prange = m_truth_info->GetParticleRange();
1087 
1088     NTruthLayers_ = 0;
1089 
1090     for (auto iter = prange.first; iter != prange.second; ++iter)
1091     {
1092         PHG4Particle *truth_particle = iter->second;
1093 
1094         // count number of layers passed through by truth particles
1095         bool layer_filled[4] = {false, false, false, false};
1096         std::map<TrkrDefs::cluskey, std::shared_ptr<TrkrCluster>> truth_clusters = truth_eval->all_truth_clusters(truth_particle);
1097         for (auto c_iter = truth_clusters.begin(); c_iter != truth_clusters.end(); c_iter++)
1098         {
1099             int layer = TrkrDefs::getLayer(c_iter->first);
1100             if (layer >= 3 && layer <= 6)
1101                 layer_filled[layer - 3] = true;
1102         }
1103         for (int i = 0; i < 4; i++)
1104             if (layer_filled[i])
1105                 NTruthLayers_++;
1106 
1107         // compute phisize and zsize of truth clusters
1108         std::set<PHG4Hit *> truth_hits = truth_eval->all_truth_hits(truth_particle);
1109         // just as with truth clustering in SvtxTruthEval, cluster together all
1110         // G4Hits in a given layer
1111         std::array<std::vector<PHG4Hit *>, 4> clusters;
1112         for (auto h_iter = truth_hits.begin(); h_iter != truth_hits.end(); ++h_iter)
1113         {
1114             PHG4Hit *hit = *h_iter;
1115             int layer = hit->get_layer();
1116             if (layer >= 3 && layer <= 6)
1117             {
1118                 clusters[layer - 3].push_back(hit);
1119             }
1120         }
1121 
1122         // since SvtxClusterEval doesn't correctly compute cluster size,
1123         // and since none of the relevant truth objects store their association
1124         // with detector elements. we have to do that association ourselves
1125         // through CylinderGeomIntt
1126         std::vector<CylinderGeomIntt *> layergeom;
1127 
1128         for (int layer = 3; layer <= 6; layer++)
1129         {
1130             layergeom.push_back(dynamic_cast<CylinderGeomIntt *>(_intt_geom_container->GetLayerGeom(layer)));
1131         }
1132 
1133         for (int i = 0; i < 4; i++)
1134         {
1135             // we break clusters up by TrkrHitset, in the same way as in
1136             // InttClusterizer
1137             std::map<TrkrDefs::hitsetkey, std::vector<PHG4Hit *>> clusters_by_hitset;
1138 
1139             for (auto &hit : clusters[i])
1140             {
1141                 double entry_point[3] = {hit->get_x(0), hit->get_y(0), hit->get_z(0)};
1142                 double exit_point[3] = {hit->get_x(1), hit->get_y(1), hit->get_z(1)};
1143 
1144                 // find ladder z and phi id
1145                 int entry_ladder_z_id, entry_ladder_phi_id;
1146                 int exit_ladder_z_id, exit_ladder_phi_id;
1147                 layergeom[i]->find_indices_from_world_location(entry_ladder_z_id, entry_ladder_phi_id, entry_point);
1148                 layergeom[i]->find_indices_from_world_location(exit_ladder_z_id, exit_ladder_phi_id, exit_point);
1149 
1150                 // fix a rounding error where you can sometimes get a phi index
1151                 // out of range, which causes a segfault when the geometry
1152                 // container can't find a surface
1153                 if (i == 0 || i == 1)
1154                 {
1155                     if (entry_ladder_phi_id == 12)
1156                         entry_ladder_phi_id = 0;
1157                     if (exit_ladder_phi_id == 12)
1158                         exit_ladder_phi_id = 0;
1159                 }
1160                 if (i == 2 || i == 3)
1161                 {
1162                     if (entry_ladder_phi_id == 16)
1163                         entry_ladder_phi_id = 0;
1164                     if (exit_ladder_phi_id == 16)
1165                         exit_ladder_phi_id = 0;
1166                 }
1167 
1168                 // get hitset key so we can retrieve surface
1169                 TrkrDefs::hitsetkey entry_hskey = InttDefs::genHitSetKey(i + 3, entry_ladder_z_id, entry_ladder_phi_id, 0);
1170                 TrkrDefs::hitsetkey exit_hskey = InttDefs::genHitSetKey(i + 3, exit_ladder_z_id, exit_ladder_phi_id, 0);
1171 
1172                 // if entry and exit ladder z id are different, then we have a
1173                 // truth G4Hit that will be split in reco (since InttClusterizer
1174                 // cannot cross TrkrHitsets) for now, just print a warning if
1175                 // this ever happens (should be very rare) and discard that hit
1176                 if (entry_hskey != exit_hskey)
1177                 {
1178                     if (Verbosity() >= VERBOSITY_MORE)
1179                     {
1180                         std::cout << "!!! WARNING !!! PHG4Hit crosses TrkrHitsets, "
1181                                      "discarding!"
1182                                   << std::endl;
1183                         std::cout << "Discarded PHG4Hit info: " << std::endl;
1184                         std::cout << "entry point: ladder (phi, z) ID = (" << entry_ladder_phi_id << ", " << entry_ladder_z_id << ")" << std::endl;
1185                         std::cout << "exit point: ladder (phi, z) ID = (" << exit_ladder_phi_id << ", " << exit_ladder_z_id << ")" << std::endl;
1186                     }
1187                     continue;
1188                 }
1189 
1190                 clusters_by_hitset[entry_hskey].push_back(hit);
1191             }
1192 
1193             for (auto chs_iter = clusters_by_hitset.begin(); chs_iter != clusters_by_hitset.end(); chs_iter++)
1194             {
1195                 TrkrDefs::hitsetkey hskey = chs_iter->first;
1196                 std::vector<PHG4Hit *> hits = chs_iter->second;
1197 
1198                 auto surface = _tgeometry->maps().getSiliconSurface(hskey);
1199 
1200                 int ladder_z_id = InttDefs::getLadderZId(hskey);
1201 
1202                 // just as in InttClusterizer, cluster size = number of unique
1203                 // strip IDs
1204                 std::set<int> phibins;
1205                 std::set<int> zbins;
1206 
1207                 std::set<TrkrDefs::cluskey> associated_reco_clusters;
1208 
1209                 for (auto &hit : hits)
1210                 {
1211                     TVector3 entry_point(hit->get_x(0), hit->get_y(0), hit->get_z(0));
1212                     TVector3 exit_point(hit->get_x(1), hit->get_y(1), hit->get_z(1));
1213 
1214                     // get local coordinates on surface
1215                     TVector3 entry_local = CylinderGeomInttHelper::get_local_from_world_coords(surface, _tgeometry, entry_point);
1216                     TVector3 exit_local = CylinderGeomInttHelper::get_local_from_world_coords(surface, _tgeometry, exit_point);
1217 
1218                     // get strip z and phi id (which we needed local coordinates
1219                     // for)
1220                     int entry_strip_phi_id, entry_strip_z_id;
1221                     int exit_strip_phi_id, exit_strip_z_id;
1222                     layergeom[i]->find_strip_index_values(ladder_z_id, entry_local[1], entry_local[2], entry_strip_phi_id, entry_strip_z_id);
1223                     layergeom[i]->find_strip_index_values(ladder_z_id, exit_local[1], exit_local[2], exit_strip_phi_id, exit_strip_z_id);
1224 
1225                     // insert one phi and z ID entry for every strip between the
1226                     // entry and exit point (this implicitly assumes that there
1227                     // will be a non-negligible amount of energy deposited along
1228                     // the whole path)
1229                     int min_z_id = std::min(entry_strip_z_id, exit_strip_z_id);
1230                     int max_z_id = std::max(entry_strip_z_id, exit_strip_z_id);
1231                     int min_phi_id = std::min(entry_strip_phi_id, exit_strip_phi_id);
1232                     int max_phi_id = std::max(entry_strip_phi_id, exit_strip_phi_id);
1233                     for (int z_id = min_z_id; z_id <= max_z_id; z_id++)
1234                     {
1235                         zbins.insert(z_id);
1236                     }
1237                     for (int phi_id = min_phi_id; phi_id <= max_phi_id; phi_id++)
1238                     {
1239                         phibins.insert(phi_id);
1240                     }
1241 
1242                     // find all reco clusters associated with this PHG4Hit
1243                     std::set<TrkrDefs::cluskey> reco_clusters_g4hit = clustereval->all_clusters_from(hit);
1244                     associated_reco_clusters.insert(reco_clusters_g4hit.begin(), reco_clusters_g4hit.end());
1245                 }
1246                 if (phibins.size() > 0)
1247                     TruthClusPhiSize_.push_back(phibins.size());
1248                 if (phibins.size() > 0 && truth_particle->get_parent_id() == 0)
1249                     PrimaryTruthClusPhiSize_.push_back(phibins.size());
1250                 if (zbins.size() > 0)
1251                     TruthClusZSize_.push_back(zbins.size());
1252                 if (zbins.size() > 0 && truth_particle->get_parent_id() == 0)
1253                     PrimaryTruthClusZSize_.push_back(zbins.size());
1254                 TruthClusNRecoClus_.push_back(associated_reco_clusters.size());
1255                 if (truth_particle->get_parent_id() == 0)
1256                     PrimaryTruthClusNRecoClus_.push_back(associated_reco_clusters.size());
1257             }
1258         }
1259     }
1260 }
1261 //____________________________________________________________________________..
1262 void dNdEtaINTT::GetPHG4Info(PHCompositeNode *topNode)
1263 {
1264     std::cout << "Get PHG4 info.: truth primary vertex" << std::endl;
1265     m_truth_info = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
1266     if (!m_truth_info)
1267     {
1268         std::cout << PHWHERE << "Error, can't find G4TruthInfo" << std::endl;
1269         exit(1);
1270     }
1271 
1272     g4hit = findNode::getClass<PHG4HitContainer>(topNode, "G4HIT_INTT");
1273     if (!g4hit)
1274     {
1275         std::cout << PHWHERE << "Error, can't find G4HIT_INTT" << std::endl;
1276         exit(1);
1277     }
1278 
1279     // Truth vertex
1280     auto vrange = m_truth_info->GetPrimaryVtxRange();
1281     int NTruthPV = 0, NTruthPV_Embeded0 = 0;
1282     for (auto iter = vrange.first; iter != vrange.second; ++iter) // process all primary vertices
1283     {
1284         const int point_id = iter->first;
1285         PHG4VtxPoint *point = iter->second;
1286         if (point)
1287         {
1288             if (m_truth_info->isEmbededVtx(point_id) == 0)
1289             {
1290                 TruthPV_trig_x_ = point->get_x();
1291                 TruthPV_trig_y_ = point->get_y();
1292                 TruthPV_trig_z_ = point->get_z();
1293                 // std::cout << "TruthVtx " << NTruthPV_Embeded0 << ", vertex: (x,y,z,t)=(" << point->get_x() << "," << point->get_y() << "," << point->get_z() << "," << point->get_t() << ")" <<
1294                 // std::endl;
1295                 NTruthPV_Embeded0++;
1296             }
1297             NTruthPV++;
1298         }
1299     }
1300     // print out the truth vertex information
1301     std::cout << "Number of truth vertices: " << NTruthPV << std::endl;
1302     std::cout << "Number of truth vertices with isEmbededVtx=0: " << NTruthPV_Embeded0 << std::endl;
1303     std::cout << "Final truth vertex: (x,y,z)=(" << TruthPV_trig_x_ << "," << TruthPV_trig_y_ << "," << TruthPV_trig_z_ << ")" << std::endl;
1304     NTruthVtx_ = NTruthPV;
1305 
1306     // PHG4Particle
1307     std::vector<int> tmpv_chargehadron;
1308     tmpv_chargehadron.clear();
1309     std::cout << "Get PHG4 info.: truth primary G4Particle" << std::endl;
1310     const auto prange = m_truth_info->GetPrimaryParticleRange();
1311     // const auto prange = m_truth_info->GetParticleRange();
1312     for (auto iter = prange.first; iter != prange.second; ++iter)
1313     {
1314         PHG4Particle *ptcl = iter->second;
1315         // particle->identify();
1316         if (ptcl)
1317         {
1318             PrimaryG4P_trackID_.push_back(ptcl->get_track_id());
1319             PrimaryG4P_PID_.push_back(ptcl->get_pid());
1320             TLorentzVector p;
1321             p.SetPxPyPzE(ptcl->get_px(), ptcl->get_py(), ptcl->get_pz(), ptcl->get_e());
1322             PrimaryG4P_E_.push_back(ptcl->get_e());
1323             PrimaryG4P_Pt_.push_back(p.Pt());
1324             PrimaryG4P_Eta_.push_back(p.Eta());
1325             PrimaryG4P_Phi_.push_back(p.Phi());
1326 
1327             TString particleclass = TString(TDatabasePDG::Instance()->GetParticle(ptcl->get_pid())->ParticleClass());
1328             bool isStable = (TDatabasePDG::Instance()->GetParticle(ptcl->get_pid())->Stable() == 1) ? true : false;
1329             double charge = TDatabasePDG::Instance()->GetParticle(ptcl->get_pid())->Charge();
1330             bool isHadron = (particleclass.Contains("Baryon") || particleclass.Contains("Meson"));
1331             bool isChargeHadron = (isStable && (charge != 0) && isHadron);
1332             if (isChargeHadron)
1333                 tmpv_chargehadron.push_back(ptcl->get_pid());
1334 
1335             PrimaryG4P_ParticleClass_.push_back(particleclass);
1336             PrimaryG4P_isStable_.push_back(isStable);
1337             PrimaryG4P_Charge_.push_back(charge);
1338             PrimaryG4P_isChargeHadron_.push_back(isChargeHadron);
1339         }
1340     }
1341     NPrimaryG4P_ = PrimaryG4P_PID_.size();
1342     NPrimaryG4P_promptChargeHadron_ = tmpv_chargehadron.size();
1343     CleanVec(tmpv_chargehadron);
1344 
1345     // PHG4Hit
1346     PHG4HitContainer::ConstRange hit_begin_end = g4hit->getHits();
1347     for (PHG4HitContainer::ConstIterator hiter = hit_begin_end.first; hiter != hit_begin_end.second; ++hiter)
1348     {
1349         PHG4Hit *hit = hiter->second;
1350         PHG4Hit_x0_.push_back(hit->get_x(0));
1351         PHG4Hit_y0_.push_back(hit->get_y(0));
1352         PHG4Hit_z0_.push_back(hit->get_z(0));
1353         PHG4Hit_x1_.push_back(hit->get_x(1));
1354         PHG4Hit_y1_.push_back(hit->get_y(1));
1355         PHG4Hit_z1_.push_back(hit->get_z(1));
1356         PHG4Hit_edep_.push_back(hit->get_edep());
1357     }
1358 }
1359 //____________________________________________________________________________..
1360 void dNdEtaINTT::GetAllPHG4Info(PHCompositeNode *topNode)
1361 {
1362     std::cout << __FILE__ << "::" << __func__ << "(): Get all PHG4 info." << std::endl;
1363 
1364     m_truth_info = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
1365     if (!m_truth_info)
1366     {
1367         std::cout << PHWHERE << "Error, can't find G4TruthInfo" << std::endl;
1368         exit(1);
1369     }
1370 
1371     const auto prange = m_truth_info->GetParticleRange();
1372     for (auto iter = prange.first; iter != prange.second; ++iter)
1373     {
1374         PHG4Particle *ptcl = iter->second;
1375         // particle->identify();
1376         if (ptcl)
1377         {
1378             G4P_PID_.push_back(ptcl->get_pid());
1379             G4P_trackID_.push_back(ptcl->get_track_id());
1380 
1381             TLorentzVector p;
1382             p.SetPxPyPzE(ptcl->get_px(), ptcl->get_py(), ptcl->get_pz(), ptcl->get_e());
1383             G4P_Pt_.push_back(p.Pt());
1384             G4P_Eta_.push_back(p.Eta());
1385             G4P_Phi_.push_back(p.Phi());
1386             G4P_E_.push_back(ptcl->get_e());
1387         }
1388     }
1389 
1390     NAllG4P_ = G4P_PID_.size();
1391 }
1392 
1393 //____________________________________________________________________________..
1394 void dNdEtaINTT::ResetVectors()
1395 {
1396     CleanVec(ClusLayer_);
1397     CleanVec(ClusX_);
1398     CleanVec(ClusY_);
1399     CleanVec(ClusZ_);
1400     CleanVec(ClusR_);
1401     CleanVec(ClusPhi_);
1402     CleanVec(ClusEta_);
1403     CleanVec(ClusAdc_);
1404     CleanVec(ClusPhiSize_);
1405     CleanVec(ClusZSize_);
1406     CleanVec(ClusLadderZId_);
1407     CleanVec(ClusLadderPhiId_);
1408     CleanVec(ClusTrkrHitSetKey_);
1409     CleanVec(ClusTimeBucketId_);
1410     CleanVec(ClusTruthCKeys_);
1411     CleanVec(ClusNG4Particles_);
1412     CleanVec(ClusNPrimaryG4Particles_);
1413     CleanVec(ClusMatchedG4P_MaxE_trackID_);
1414     CleanVec(ClusMatchedG4P_MaxE_Pt_);
1415     CleanVec(ClusMatchedG4P_MaxE_Eta_);
1416     CleanVec(ClusMatchedG4P_MaxE_Phi_);
1417     CleanVec(ClusMatchedG4P_MaxClusE_trackID_);
1418     CleanVec(ClusMatchedG4P_MaxClusE_ancestorTrackID_);
1419     CleanVec(ClusMatchedG4P_MaxClusE_Pt_);
1420     CleanVec(ClusMatchedG4P_MaxClusE_Eta_);
1421     CleanVec(ClusMatchedG4P_MaxClusE_Phi_);
1422     CleanVec(TruthClusPhiSize_);
1423     CleanVec(TruthClusZSize_);
1424     CleanVec(TruthClusNRecoClus_);
1425     CleanVec(PrimaryTruthClusPhiSize_);
1426     CleanVec(PrimaryTruthClusZSize_);
1427     CleanVec(PrimaryTruthClusNRecoClus_);
1428     CleanVec(InttRawHit_bco_);
1429     CleanVec(InttRawHit_packetid_);
1430     CleanVec(InttRawHit_word_);
1431     CleanVec(InttRawHit_fee_);
1432     CleanVec(InttRawHit_channel_id_);
1433     CleanVec(InttRawHit_chip_id_);
1434     CleanVec(InttRawHit_adc_);
1435     CleanVec(InttRawHit_FPHX_BCO_);
1436     CleanVec(InttRawHit_full_FPHX_);
1437     CleanVec(InttRawHit_full_ROC_);
1438     CleanVec(InttRawHit_amplitude_);
1439     CleanVec(TrkrHitRow_);
1440     CleanVec(TrkrHitColumn_);
1441     CleanVec(TrkrHitLadderZId_);
1442     CleanVec(TrkrHitLadderPhiId_);
1443     CleanVec(TrkrHitTimeBucketId_);
1444     CleanVec(TrkrHitLayer_);
1445     CleanVec(TrkrHitADC_);
1446     CleanVec(TrkrHitX_);
1447     CleanVec(TrkrHitY_);
1448     CleanVec(TrkrHitZ_);
1449     CleanVec(TrkrHit_truthHit_x0_);
1450     CleanVec(TrkrHit_truthHit_y0_);
1451     CleanVec(TrkrHit_truthHit_z0_);
1452     CleanVec(TrkrHit_truthHit_x1_);
1453     CleanVec(TrkrHit_truthHit_y1_);
1454     CleanVec(TrkrHit_truthHit_z1_);
1455     CleanVec(HepMCFSPrtl_Pt_);
1456     CleanVec(HepMCFSPrtl_Eta_);
1457     CleanVec(HepMCFSPrtl_Phi_);
1458     CleanVec(HepMCFSPrtl_E_);
1459     CleanVec(HepMCFSPrtl_PID_);
1460     CleanVec(HepMCFSPrtl_prodx_);
1461     CleanVec(HepMCFSPrtl_prody_);
1462     CleanVec(HepMCFSPrtl_prodz_);
1463     CleanVec(PrimaryG4P_Pt_);
1464     CleanVec(PrimaryG4P_Eta_);
1465     CleanVec(PrimaryG4P_Phi_);
1466     CleanVec(PrimaryG4P_E_);
1467     CleanVec(PrimaryG4P_PID_);
1468     CleanVec(PrimaryG4P_trackID_);
1469     CleanVec(PrimaryG4P_ParticleClass_);
1470     CleanVec(PrimaryG4P_isStable_);
1471     CleanVec(PrimaryG4P_Charge_);
1472     CleanVec(PrimaryG4P_isChargeHadron_);
1473     CleanVec(PHG4Hit_x0_);
1474     CleanVec(PHG4Hit_y0_);
1475     CleanVec(PHG4Hit_z0_);
1476     CleanVec(PHG4Hit_x1_);
1477     CleanVec(PHG4Hit_y1_);
1478     CleanVec(PHG4Hit_z1_);
1479     CleanVec(PHG4Hit_edep_);
1480     CleanVec(G4P_Pt_);
1481     CleanVec(G4P_Eta_);
1482     CleanVec(G4P_Phi_);
1483     CleanVec(G4P_E_);
1484     CleanVec(G4P_PID_);
1485     CleanVec(G4P_trackID_);
1486     CleanVec(firedTriggers_);
1487     CleanVec(firedTriggers_name_);
1488     CleanVec(firedTriggers_checkraw_);
1489     CleanVec(firedTriggers_prescale_);
1490     CleanVec(firedTriggers_scalers_);
1491     CleanVec(firedTriggers_livescalers_);
1492     CleanVec(firedTriggers_rawscalers_);
1493 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!