ChargedHadronMultiplicity-SmallSystem/module/VertexCompare.cc

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 VertexCompare.h.
0007 //
0008 // VertexCompare(const std::string &name = "VertexCompare")
0009 // everything is keyed to VertexCompare, 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 // VertexCompare::~VertexCompare()
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 node tree
0016 //
0017 // int VertexCompare::Init(PHCompositeNode *topNode)
0018 // This method is called when the module is registered with the Fun4AllServer. You
0019 // can create historgrams here or put objects on the node tree but be aware that
0020 // modules which haven't been registered yet did not put antyhing on the node tree
0021 //
0022 // int VertexCompare::InitRun(PHCompositeNode *topNode)
0023 // This method is called when the first event is read (or generated). At
0024 // this point the run number is known (which is mainly interesting for raw data
0025 // processing). Also all objects are on the node tree in case your module's action
0026 // depends on what else is around. Last chance to put nodes under the DST Node
0027 // We mix events during readback if branches are added after the first event
0028 //
0029 // int VertexCompare::process_event(PHCompositeNode *topNode)
0030 // called for every event. Return codes trigger actions, you find them in
0031 // $OFFLINE_MAIN/include/fun4all/Fun4AllReturnCodes.h
0032 //   everything is good:
0033 //     return Fun4AllReturnCodes::EVENT_OK
0034 //   abort event reconstruction, clear everything and process next event:
0035 //     return Fun4AllReturnCodes::ABORT_EVENT;
0036 //   proceed but do not save this event in output (needs output manager setting):
0037 //     return Fun4AllReturnCodes::DISCARD_EVENT;
0038 //   abort processing:
0039 //     return Fun4AllReturnCodes::ABORT_RUN
0040 // all other integers will lead to an error and abort of processing
0041 //
0042 // int VertexCompare::ResetEvent(PHCompositeNode *topNode)
0043 // If you have internal data structures (arrays, stl containers) which needs clearing
0044 // after each event, this is the place to do that. The nodes under the DST node are cleared
0045 // by the framework
0046 //
0047 // int VertexCompare::EndRun(const int runnumber)
0048 // This method is called at the end of a run when an event from a new run is
0049 // encountered. Useful when analyzing multiple runs (raw data). Also called at
0050 // the end of processing (before the End() method)
0051 //
0052 // int VertexCompare::End(PHCompositeNode *topNode)
0053 // This is called at the end of processing. It needs to be called by the macro
0054 // by Fun4AllServer::End(), so do not forget this in your macro
0055 //
0056 // int VertexCompare::Reset(PHCompositeNode *topNode)
0057 // not really used - it is called before the dtor is called
0058 //
0059 // void VertexCompare::Print(const std::string &what) const
0060 // Called from the command line - useful to print information when you need it
0061 //
0062 //____________________________________________________________________________..
0063
0064 #include "VertexCompare.h"
0065
0066 #include <fun4all/Fun4AllReturnCodes.h>
0067
0068 #include <phool/PHCompositeNode.h>
0069 #include <phool/getClass.h>
0070
0071 #include <trackbase/ActsGeometry.h>
0072 #include <trackbase/TrkrCluster.h>
0073 #include <trackbase/TrkrClusterContainer.h>
0074 #include <trackbase/TrkrClusterHitAssoc.h>
0075 #include <trackbase_historic/SvtxTrack.h>
0076 #include <trackbase_historic/SvtxTrackMap.h>
0077 #include <trackbase_historic/TrackAnalysisUtils.h>
0078 #include <trackbase_historic/TrackSeed.h>
0079 #include <trackbase_historic/TrackSeedContainer.h>
0080 #include <trackbase_historic/TrackSeedHelper.h>
0081
0082 #include <mvtx/SegmentationAlpide.h>
0083 #include <mvtx/MvtxPixelDefs.h>
0084
0085 #include <trackbase/InttDefs.h>
0086 #include <trackbase/MvtxDefs.h>
0087 #include <trackbase/TrkrDefs.h>
0088
0089 #include <ffarawobjects/Gl1Packet.h>
0090 #include <ffarawobjects/Gl1RawHit.h>
0091
0092 #include <globalvertex/GlobalVertex.h>
0093 #include <globalvertex/GlobalVertexMap.h>
0094 #include <globalvertex/MbdVertex.h>
0095 #include <globalvertex/MbdVertexMap.h>
0096 #include <globalvertex/SvtxVertex.h>
0097 #include <globalvertex/SvtxVertexMap.h>
0098
0099 #include <mbd/MbdOut.h>
0100 #include <mbd/MbdPmtContainer.h>
0101 #include <mbd/MbdPmtHit.h>
0102
0103 #include <calotrigger/MinimumBiasInfo.h>
0104 #include <centrality/CentralityInfo.h>
0105
0106 #include "g4eval/SvtxClusterEval.h"
0107 #include <g4eval/SvtxEvalStack.h>
0108 #include <g4eval/SvtxHitEval.h>
0109 #include <g4eval/SvtxTruthEval.h>
0110 #include <g4main/PHG4Hit.h>
0111 #include <g4main/PHG4HitContainer.h>
0112 #include <g4main/PHG4Particle.h>
0113 #include <g4main/PHG4MCProcessDefs.h>
0114 #include <g4main/PHG4TruthInfoContainer.h>
0115 #include <g4main/PHG4VtxPoint.h>
0116
0117 #include <algorithm>
0118 #include <map>
0119
0120 #include <TDatabasePDG.h>
0121 #include <TParticlePDG.h>
0122
0123 namespace
0124 {
0125 template <class Container> void Clean(Container &c) { Container().swap(c); }
0126 } // namespace
0127
0128 GlobalVertex::VTXTYPE trkType = GlobalVertex::SVTX;
0129 GlobalVertex::VTXTYPE mbdType = GlobalVertex::MBD;
0130 //____________________________________________________________________________..
0131 VertexCompare::VertexCompare(const std::string &name)
0132     : SubsysReco(name)
0133 {
0134     std::cout << "VertexCompare::VertexCompare(const std::string &name) Calling ctor" << std::endl;
0135 }
0136
0137 //____________________________________________________________________________..
0138 VertexCompare::~VertexCompare() { std::cout << "VertexCompare::~VertexCompare() Calling dtor" << std::endl; }
0139
0140 //____________________________________________________________________________..
0141 int VertexCompare::Init(PHCompositeNode *topNode)
0142 {
0143     outFile = new TFile(outFileName.c_str(), "RECREATE");
0144     outTree = new TTree("VTX", "VTX");
0145     // outTree->OptimizeBaskets();
0146     // outTree->SetAutoSave(-5e6);
0147
0148     outTree->Branch("counter", &counter, "counter/I");
0149     outTree->Branch("is_min_bias", &is_min_bias);
0150     outTree->Branch("firedTriggers", &firedTriggers);
0151     outTree->Branch("gl1bco", &gl1bco);
0152     outTree->Branch("gl1BunchCrossing", &gl1BunchCrossing);
0153     outTree->Branch("bcotr", &bcotr);
0154     outTree->Branch("centrality_mbd", &centrality_mbd_);
0155     outTree->Branch("n_MBDVertex", &n_MBDVertex, "n_MBDVertex/i");
0156     outTree->Branch("mbdVertex", &mbdVertex);
0157     outTree->Branch("mbdVertexId", &mbdVertexId);
0158     outTree->Branch("mbdVertexCrossing", &mbdVertexCrossing);
0159     outTree->Branch("MBD_charge_sum", &MBD_charge_sum);
0160     outTree->Branch("nSvtxVertices", &nSvtxVertices);
0161     outTree->Branch("trackerVertexId", &trackerVertexId);
0162     outTree->Branch("trackerVertexX", &trackerVertexX);
0163     outTree->Branch("trackerVertexY", &trackerVertexY);
0164     outTree->Branch("trackerVertexZ", &trackerVertexZ);
0165     outTree->Branch("trackerVertexChisq", &trackerVertexChisq);
0166     outTree->Branch("trackerVertexNdof", &trackerVertexNdof);
0167     outTree->Branch("trackerVertexNTracks", &trackerVertexNTracks);
0168     outTree->Branch("trackerVertexCrossing", &trackerVertexCrossing);
0169     outTree->Branch("trackerVertexTrackIDs", &trackerVertexTrackIDs);
0170     outTree->Branch("nTracks", &nTracks, "nTracks/i");
0171     // outTree->Branch("n_TRKVertex", &n_TRKVertex, "n_TRKVertex/i");
0172     outTree->Branch("hasMBD", &hasMBD, "hasMBD/O");
0173     outTree->Branch("hasTRK", &hasTRK, "hasTRK/O");
0174
0175     // silicon seed information
0176     outTree->Branch("nTotalSilSeeds", &nTotalSilSeeds);
0177     outTree->Branch("nSilSeedsValidCrossing", &nSilSeedsValidCrossing);
0178     outTree->Branch("silseed_id", &silseed_id);
0179     outTree->Branch("silseed_assocVtxId", &silseed_assocVtxId);
0180     outTree->Branch("silseed_x", &silseed_x);
0181     outTree->Branch("silseed_y", &silseed_y);
0182     outTree->Branch("silseed_z", &silseed_z);
0183     outTree->Branch("silseed_pt", &silseed_pt);
0184     outTree->Branch("silseed_eta", &silseed_eta);
0185     outTree->Branch("silseed_phi", &silseed_phi);
0186     outTree->Branch("silseed_eta_vtx", &silseed_eta_vtx);
0187     outTree->Branch("silseed_phi_vtx", &silseed_phi_vtx);
0188     outTree->Branch("silseed_crossing", &silseed_crossing);
0189     outTree->Branch("silseed_charge", &silseed_charge);
0190     outTree->Branch("silseed_nMvtx", &silseed_nMvtx);
0191     outTree->Branch("silseed_nIntt", &silseed_nIntt);
0192     outTree->Branch("silseed_clusterKeys", &silseed_clusterKeys);
0193     outTree->Branch("silseed_cluster_layer", &silseed_cluster_layer);
0194     outTree->Branch("silseed_cluster_globalX", &silseed_cluster_globalX);
0195     outTree->Branch("silseed_cluster_globalY", &silseed_cluster_globalY);
0196     outTree->Branch("silseed_cluster_globalZ", &silseed_cluster_globalZ);
0197     outTree->Branch("silseed_cluster_phi", &silseed_cluster_phi);
0198     outTree->Branch("silseed_cluster_eta", &silseed_cluster_eta);
0199     outTree->Branch("silseed_cluster_r", &silseed_cluster_r);
0200     outTree->Branch("silseed_cluster_phiSize", &silseed_cluster_phiSize);
0201     outTree->Branch("silseed_cluster_zSize", &silseed_cluster_zSize);
0202     outTree->Branch("silseed_cluster_strobeID", &silseed_cluster_strobeID);
0203     outTree->Branch("silseed_cluster_timeBucketID", &silseed_cluster_timeBucketID); // only for INTT
0204
0205     // (intt) cluster information
0206     outTree->Branch("clusterKey", &clusterKey);
0207     outTree->Branch("cluster_layer", &cluster_layer);
0208     outTree->Branch("cluster_chip", &cluster_chip);   // for mvtx chip id
0209     outTree->Branch("cluster_stave", &cluster_stave); // for mvtx stave id
0210     outTree->Branch("cluster_globalX", &cluster_globalX);
0211     outTree->Branch("cluster_globalY", &cluster_globalY);
0212     outTree->Branch("cluster_globalZ", &cluster_globalZ);
0213     outTree->Branch("cluster_phi", &cluster_phi);
0214     outTree->Branch("cluster_eta", &cluster_eta);
0215     outTree->Branch("cluster_r", &cluster_r);
0216     outTree->Branch("cluster_phiSize", &cluster_phiSize);
0217     outTree->Branch("cluster_zSize", &cluster_zSize);
0218     outTree->Branch("cluster_adc", &cluster_adc);
0219     outTree->Branch("cluster_timeBucketID", &cluster_timeBucketID);
0220     outTree->Branch("cluster_ladderZId", &cluster_ladderZId);     // for intt ladder z id
0221     outTree->Branch("cluster_ladderPhiId", &cluster_ladderPhiId); // for intt ladder phi id
0222     outTree->Branch("cluster_LocalX", &cluster_LocalX);
0223     outTree->Branch("cluster_LocalY", &cluster_LocalY);
0224     // cluster truth matching by max_truth_particle_by_energy
0225     outTree->Branch("cluster_matchedG4P_trackID", &cluster_matchedG4P_trackID);
0226     outTree->Branch("cluster_matchedG4P_PID", &cluster_matchedG4P_PID);
0227     outTree->Branch("cluster_matchedG4P_E", &cluster_matchedG4P_E);
0228     outTree->Branch("cluster_matchedG4P_pT", &cluster_matchedG4P_pT);
0229     outTree->Branch("cluster_matchedG4P_eta", &cluster_matchedG4P_eta);
0230     outTree->Branch("cluster_matchedG4P_phi", &cluster_matchedG4P_phi);
0231
0232     outTree->Branch("mvtx_seedcluster_key", &mvtx_seedcluster_key);
0233     outTree->Branch("mvtx_seedcluster_layer", &mvtx_seedcluster_layer);
0234     outTree->Branch("mvtx_seedcluster_chip", &mvtx_seedcluster_chip);
0235     outTree->Branch("mvtx_seedcluster_stave", &mvtx_seedcluster_stave);
0236     outTree->Branch("mvtx_seedcluster_globalX", &mvtx_seedcluster_globalX);
0237     outTree->Branch("mvtx_seedcluster_globalY", &mvtx_seedcluster_globalY);
0238     outTree->Branch("mvtx_seedcluster_globalZ", &mvtx_seedcluster_globalZ);
0239     outTree->Branch("mvtx_seedcluster_phi", &mvtx_seedcluster_phi);
0240     outTree->Branch("mvtx_seedcluster_eta", &mvtx_seedcluster_eta);
0241     outTree->Branch("mvtx_seedcluster_r", &mvtx_seedcluster_r);
0242     outTree->Branch("mvtx_seedcluster_phiSize", &mvtx_seedcluster_phiSize);
0243     outTree->Branch("mvtx_seedcluster_zSize", &mvtx_seedcluster_zSize);
0244     outTree->Branch("mvtx_seedcluster_strobeID", &mvtx_seedcluster_strobeID);
0245     outTree->Branch("mvtx_seedcluster_matchedcrossing", &mvtx_seedcluster_matchedcrossing);
0246     outTree->Branch("mvtx_seedcluster_hitX", &mvtx_seedcluster_hitX);
0247     outTree->Branch("mvtx_seedcluster_hitY", &mvtx_seedcluster_hitY);
0248     outTree->Branch("mvtx_seedcluster_hitZ", &mvtx_seedcluster_hitZ);
0249     outTree->Branch("mvtx_seedcluster_hitrow", &mvtx_seedcluster_hitrow);
0250     outTree->Branch("mvtx_seedcluster_hitcol", &mvtx_seedcluster_hitcol);
0251
0252     if (isSimulation)
0253     {
0254         outTree->Branch("nTruthVertex", &nTruthVertex);
0255         outTree->Branch("TruthVertexX", &TruthVertexX);
0256         outTree->Branch("TruthVertexY", &TruthVertexY);
0257         outTree->Branch("TruthVertexZ", &TruthVertexZ);
0258
0259         outTree->Branch("silseed_ngmvtx", &silseed_ngmvtx);
0260         outTree->Branch("silseed_ngintt", &silseed_ngintt);
0261         outTree->Branch("silseed_cluster_gcluster_key", &silseed_cluster_gcluster_key);
0262         outTree->Branch("silseed_cluster_gcluster_layer", &silseed_cluster_gcluster_layer);
0263         outTree->Branch("silseed_cluster_gcluster_X", &silseed_cluster_gcluster_X);
0264         outTree->Branch("silseed_cluster_gcluster_Y", &silseed_cluster_gcluster_Y);
0265         outTree->Branch("silseed_cluster_gcluster_Z", &silseed_cluster_gcluster_Z);
0266         outTree->Branch("silseed_cluster_gcluster_r", &silseed_cluster_gcluster_r);
0267         outTree->Branch("silseed_cluster_gcluster_phi", &silseed_cluster_gcluster_phi);
0268         outTree->Branch("silseed_cluster_gcluster_eta", &silseed_cluster_gcluster_eta);
0269         outTree->Branch("silseed_cluster_gcluster_edep", &silseed_cluster_gcluster_edep);
0270         outTree->Branch("silseed_cluster_gcluster_adc", &silseed_cluster_gcluster_adc);
0271         outTree->Branch("silseed_cluster_gcluster_phiSize", &silseed_cluster_gcluster_phiSize);
0272         outTree->Branch("silseed_cluster_gcluster_zSize", &silseed_cluster_gcluster_zSize);
0273
0274         outTree->Branch("mvtx_seedcluster_matchedG4P_trackID", &mvtx_seedcluster_matchedG4P_trackID);
0275         outTree->Branch("mvtx_seedcluster_matchedG4P_PID", &mvtx_seedcluster_matchedG4P_PID);
0276         outTree->Branch("mvtx_seedcluster_matchedG4P_E", &mvtx_seedcluster_matchedG4P_E);
0277         outTree->Branch("mvtx_seedcluster_matchedG4P_pT", &mvtx_seedcluster_matchedG4P_pT);
0278         outTree->Branch("mvtx_seedcluster_matchedG4P_eta", &mvtx_seedcluster_matchedG4P_eta);
0279         outTree->Branch("mvtx_seedcluster_matchedG4P_phi", &mvtx_seedcluster_matchedG4P_phi);
0280         outTree->Branch("mvtx_seedcluster_matchedG4P_ancestor_trackID", &mvtx_seedcluster_matchedG4P_ancestor_trackID);
0281         outTree->Branch("mvtx_seedcluster_matchedG4P_ancestor_PID", &mvtx_seedcluster_matchedG4P_ancestor_PID);
0282
0283         outTree->Branch("N_PrimaryPHG4Ptcl", &N_PrimaryPHG4Ptcl);
0284         outTree->Branch("N_sPHENIXPrimary", &N_sPHENIXPrimary);
0285         outTree->Branch("N_AllPHG4Ptcl", &N_AllPHG4Ptcl);
0286
0287         outTree->Branch("PrimaryPHG4Ptcl_pT", &PrimaryPHG4Ptcl_pT);
0288         outTree->Branch("PrimaryPHG4Ptcl_eta", &PrimaryPHG4Ptcl_eta);
0289         outTree->Branch("PrimaryPHG4Ptcl_phi", &PrimaryPHG4Ptcl_phi);
0290         outTree->Branch("PrimaryPHG4Ptcl_E", &PrimaryPHG4Ptcl_E);
0291         outTree->Branch("PrimaryPHG4Ptcl_PID", &PrimaryPHG4Ptcl_PID);
0292         outTree->Branch("PrimaryPHG4Ptcl_trackID", &PrimaryPHG4Ptcl_trackID);
0293         outTree->Branch("PrimaryPHG4Ptcl_ParticleClass", &PrimaryPHG4Ptcl_ParticleClass);
0294         outTree->Branch("PrimaryPHG4Ptcl_isStable", &PrimaryPHG4Ptcl_isStable);
0295         outTree->Branch("PrimaryPHG4Ptcl_charge", &PrimaryPHG4Ptcl_charge);
0296         outTree->Branch("PrimaryPHG4Ptcl_isChargedHadron", &PrimaryPHG4Ptcl_isChargedHadron);
0297         outTree->Branch("PrimaryPHG4Ptcl_ancestor_trackID", &PrimaryPHG4Ptcl_ancestor_trackID);
0298         outTree->Branch("PrimaryPHG4Ptcl_ancestor_PID", &PrimaryPHG4Ptcl_ancestor_PID);
0299         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_X", &PrimaryPHG4Ptcl_truthcluster_X);
0300         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_Y", &PrimaryPHG4Ptcl_truthcluster_Y);
0301         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_Z", &PrimaryPHG4Ptcl_truthcluster_Z);
0302         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_edep", &PrimaryPHG4Ptcl_truthcluster_edep);
0303         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_adc", &PrimaryPHG4Ptcl_truthcluster_adc);
0304         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_r", &PrimaryPHG4Ptcl_truthcluster_r);
0305         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_phi", &PrimaryPHG4Ptcl_truthcluster_phi);
0306         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_eta", &PrimaryPHG4Ptcl_truthcluster_eta);
0307         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_phisize", &PrimaryPHG4Ptcl_truthcluster_phisize);
0308         outTree->Branch("PrimaryPHG4Ptcl_truthcluster_zsize", &PrimaryPHG4Ptcl_truthcluster_zsize);
0309         outTree->Branch("PrimaryPHG4Ptcl_recocluster_globalX", &PrimaryPHG4Ptcl_recocluster_globalX);
0310         outTree->Branch("PrimaryPHG4Ptcl_recocluster_globalY", &PrimaryPHG4Ptcl_recocluster_globalY);
0311         outTree->Branch("PrimaryPHG4Ptcl_recocluster_globalZ", &PrimaryPHG4Ptcl_recocluster_globalZ);
0312         outTree->Branch("PrimaryPHG4Ptcl_recocluster_r", &PrimaryPHG4Ptcl_recocluster_r);
0313         outTree->Branch("PrimaryPHG4Ptcl_recocluster_phi", &PrimaryPHG4Ptcl_recocluster_phi);
0314         outTree->Branch("PrimaryPHG4Ptcl_recocluster_eta", &PrimaryPHG4Ptcl_recocluster_eta);
0315         outTree->Branch("PrimaryPHG4Ptcl_recocluster_phisize", &PrimaryPHG4Ptcl_recocluster_phisize);
0316         outTree->Branch("PrimaryPHG4Ptcl_recocluster_zsize", &PrimaryPHG4Ptcl_recocluster_zsize);
0317         outTree->Branch("PrimaryPHG4Ptcl_recocluster_adc", &PrimaryPHG4Ptcl_recocluster_adc);
0318
0319         outTree->Branch("sPHENIXPrimary_pT", &sPHENIXPrimary_pT);
0320         outTree->Branch("sPHENIXPrimary_eta", &sPHENIXPrimary_eta);
0321         outTree->Branch("sPHENIXPrimary_phi", &sPHENIXPrimary_phi);
0322         outTree->Branch("sPHENIXPrimary_E", &sPHENIXPrimary_E);
0323         outTree->Branch("sPHENIXPrimary_PID", &sPHENIXPrimary_PID);
0324         outTree->Branch("sPHENIXPrimary_trackID", &sPHENIXPrimary_trackID);
0325         outTree->Branch("sPHENIXPrimary_ParticleClass", &sPHENIXPrimary_ParticleClass);
0326         outTree->Branch("sPHENIXPrimary_isStable", &sPHENIXPrimary_isStable);
0327         outTree->Branch("sPHENIXPrimary_charge", &sPHENIXPrimary_charge);
0328         outTree->Branch("sPHENIXPrimary_isChargedHadron", &sPHENIXPrimary_isChargedHadron);
0329         outTree->Branch("sPHENIXPrimary_ancestor_trackID", &sPHENIXPrimary_ancestor_trackID);
0330         outTree->Branch("sPHENIXPrimary_ancestor_PID", &sPHENIXPrimary_ancestor_PID);
0331         outTree->Branch("sPHENIXPrimary_truthcluster_X", &sPHENIXPrimary_truthcluster_X);
0332         outTree->Branch("sPHENIXPrimary_truthcluster_Y", &sPHENIXPrimary_truthcluster_Y);
0333         outTree->Branch("sPHENIXPrimary_truthcluster_Z", &sPHENIXPrimary_truthcluster_Z);
0334         outTree->Branch("sPHENIXPrimary_truthcluster_edep", &sPHENIXPrimary_truthcluster_edep);
0335         outTree->Branch("sPHENIXPrimary_truthcluster_adc", &sPHENIXPrimary_truthcluster_adc);
0336         outTree->Branch("sPHENIXPrimary_truthcluster_r", &sPHENIXPrimary_truthcluster_r);
0337         outTree->Branch("sPHENIXPrimary_truthcluster_phi", &sPHENIXPrimary_truthcluster_phi);
0338         outTree->Branch("sPHENIXPrimary_truthcluster_eta", &sPHENIXPrimary_truthcluster_eta);
0339         outTree->Branch("sPHENIXPrimary_truthcluster_phisize", &sPHENIXPrimary_truthcluster_phisize);
0340         outTree->Branch("sPHENIXPrimary_truthcluster_zsize", &sPHENIXPrimary_truthcluster_zsize);
0341         outTree->Branch("sPHENIXPrimary_recocluster_globalX", &sPHENIXPrimary_recocluster_globalX);
0342         outTree->Branch("sPHENIXPrimary_recocluster_globalY", &sPHENIXPrimary_recocluster_globalY);
0343         outTree->Branch("sPHENIXPrimary_recocluster_globalZ", &sPHENIXPrimary_recocluster_globalZ);
0344         outTree->Branch("sPHENIXPrimary_recocluster_r", &sPHENIXPrimary_recocluster_r);
0345         outTree->Branch("sPHENIXPrimary_recocluster_phi", &sPHENIXPrimary_recocluster_phi);
0346         outTree->Branch("sPHENIXPrimary_recocluster_eta", &sPHENIXPrimary_recocluster_eta);
0347         outTree->Branch("sPHENIXPrimary_recocluster_phisize", &sPHENIXPrimary_recocluster_phisize);
0348         outTree->Branch("sPHENIXPrimary_recocluster_zsize", &sPHENIXPrimary_recocluster_zsize);
0349         outTree->Branch("sPHENIXPrimary_recocluster_adc", &sPHENIXPrimary_recocluster_adc);
0350
0351         outTree->Branch("AllPHG4Ptcl_pT", &AllPHG4Ptcl_pT);
0352         outTree->Branch("AllPHG4Ptcl_eta", &AllPHG4Ptcl_eta);
0353         outTree->Branch("AllPHG4Ptcl_phi", &AllPHG4Ptcl_phi);
0354         outTree->Branch("AllPHG4Ptcl_E", &AllPHG4Ptcl_E);
0355         outTree->Branch("AllPHG4Ptcl_PID", &AllPHG4Ptcl_PID);
0356         outTree->Branch("AllPHG4Ptcl_trackID", &AllPHG4Ptcl_trackID);
0357         outTree->Branch("AllPHG4Ptcl_ancestor_trackID", &AllPHG4Ptcl_ancestor_trackID);
0358         outTree->Branch("AllPHG4Ptcl_ancestor_PID", &AllPHG4Ptcl_ancestor_PID);
0359         outTree->Branch("AllPHG4Ptcl_truthcluster_X", &AllPHG4Ptcl_truthcluster_X);
0360         outTree->Branch("AllPHG4Ptcl_truthcluster_Y", &AllPHG4Ptcl_truthcluster_Y);
0361         outTree->Branch("AllPHG4Ptcl_truthcluster_Z", &AllPHG4Ptcl_truthcluster_Z);
0362         outTree->Branch("AllPHG4Ptcl_truthcluster_edep", &AllPHG4Ptcl_truthcluster_edep);
0363         outTree->Branch("AllPHG4Ptcl_truthcluster_adc", &AllPHG4Ptcl_truthcluster_adc);
0364         outTree->Branch("AllPHG4Ptcl_truthcluster_r", &AllPHG4Ptcl_truthcluster_r);
0365         outTree->Branch("AllPHG4Ptcl_truthcluster_phi", &AllPHG4Ptcl_truthcluster_phi);
0366         outTree->Branch("AllPHG4Ptcl_truthcluster_eta", &AllPHG4Ptcl_truthcluster_eta);
0367         outTree->Branch("AllPHG4Ptcl_truthcluster_phisize", &AllPHG4Ptcl_truthcluster_phisize);
0368         outTree->Branch("AllPHG4Ptcl_truthcluster_zsize", &AllPHG4Ptcl_truthcluster_zsize);
0369         outTree->Branch("AllPHG4Ptcl_recocluster_globalX", &AllPHG4Ptcl_recocluster_globalX);
0370         outTree->Branch("AllPHG4Ptcl_recocluster_globalY", &AllPHG4Ptcl_recocluster_globalY);
0371         outTree->Branch("AllPHG4Ptcl_recocluster_globalZ", &AllPHG4Ptcl_recocluster_globalZ);
0372         outTree->Branch("AllPHG4Ptcl_recocluster_r", &AllPHG4Ptcl_recocluster_r);
0373         outTree->Branch("AllPHG4Ptcl_recocluster_phi", &AllPHG4Ptcl_recocluster_phi);
0374         outTree->Branch("AllPHG4Ptcl_recocluster_eta", &AllPHG4Ptcl_recocluster_eta);
0375         outTree->Branch("AllPHG4Ptcl_recocluster_phisize", &AllPHG4Ptcl_recocluster_phisize);
0376         outTree->Branch("AllPHG4Ptcl_recocluster_zsize", &AllPHG4Ptcl_recocluster_zsize);
0377         outTree->Branch("AllPHG4Ptcl_recocluster_adc", &AllPHG4Ptcl_recocluster_adc);
0378     }
0379
0380     return Fun4AllReturnCodes::EVENT_OK;
0381 }
0382
0383 //____________________________________________________________________________..
0384 int VertexCompare::InitRun(PHCompositeNode *topNode) { return Fun4AllReturnCodes::EVENT_OK; }
0385
0386 //____________________________________________________________________________..
0387 int VertexCompare::process_event(PHCompositeNode *topNode)
0388 {
0389     std::cout << "VertexCompare::process_event - Processing event " << counter << std::endl;
0390
0391     PHNodeIterator dstiter(topNode);
0392     PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(dstiter.findFirst("PHCompositeNode", "DST"));
0393
0394     MbdVertexMap *m_dst_mbdvertexmap = findNode::getClass<MbdVertexMap>(topNode, "MbdVertexMap");
0395     SvtxVertexMap *m_dst_vertexmap = findNode::getClass<SvtxVertexMap>(topNode, "SvtxVertexMap");
0396
0397     auto globalvertexmap = findNode::getClass<GlobalVertexMap>(topNode, "GlobalVertexMap");
0398
0399     clustermap = findNode::getClass<TrkrClusterContainer>(dstNode, clusterContainerName);
0400     clusterhitassoc = findNode::getClass<TrkrClusterHitAssoc>(dstNode, clusterHitAssocName);
0401     geometry = findNode::getClass<ActsGeometry>(topNode, geometryNodeName);
0402     silseedmap = findNode::getClass<TrackSeedContainer>(topNode, seedContainerName);
0403     gl1PacketInfo = findNode::getClass<Gl1Packet>(topNode, gl1NodeName);
0404     m_mbdout = findNode::getClass<MbdOut>(topNode, mbdOutNodeName);
0405     minimumbiasinfo = findNode::getClass<MinimumBiasInfo>(topNode, "MinimumBiasInfo");
0406     m_CentInfo = findNode::getClass<CentralityInfo>(topNode, "CentralityInfo");
0407
0408     if (!clustermap)
0409     {
0410         std::cout << "SiliconSeedAnalyzer::process_event - [ERROR] - can't find cluster map node " << clusterContainerName << std::endl;
0411         // return Fun4AllReturnCodes::ABORTEVENT;
0412     }
0413     if (!clusterhitassoc)
0414     {
0415         std::cout << "SiliconSeedAnalyzer::process_event - [ERROR] - can't find cluster hit association node " << clusterHitAssocName << std::endl;
0416         // return Fun4AllReturnCodes::ABORTEVENT;
0417     }
0418     if (!geometry)
0419     {
0420         std::cout << "SiliconSeedAnalyzer::process_event - [ERROR] - can't find ActsGeometry node " << geometryNodeName << std::endl;
0421         // return Fun4AllReturnCodes::ABORTEVENT;
0422     }
0423     if (!silseedmap)
0424     {
0425         std::cout << "SiliconSeedAnalyzer::process_event - [ERROR] - can't find silicon seed map node " << seedContainerName << std::endl;
0426         // return Fun4AllReturnCodes::ABORTEVENT;
0427     }
0428     if (!gl1PacketInfo)
0429     {
0430         std::cout << "SiliconSeedAnalyzer::process_event - [WARNING] - can't find GL1 node " << gl1NodeName << std::endl;
0431         // return Fun4AllReturnCodes::EVENT_OK;
0432     }
0433     else
0434     {
0435         gl1BunchCrossing = gl1PacketInfo->getBunchNumber();
0436         uint64_t triggervec = gl1PacketInfo->getScaledVector();
0437         gl1bco = gl1PacketInfo->getBCO();
0438         auto lbshift = gl1bco << 24U;
0439         bcotr = lbshift >> 24U;
0440         for (int i = 0; i < 64; ++i)
0441         {
0442             bool trig_decision = ((triggervec & 0x1U) == 0x1U);
0443             if (trig_decision)
0444             {
0445                 firedTriggers.push_back(i);
0446             }
0447             triggervec = (triggervec >> 1U) & 0xffffffffU;
0448         }
0449     }
0450     if (!m_mbdout)
0451     {
0452         std::cout << "SiliconSeedAnalyzer::process_event - [WARNING] - can't find MbdOut node " << mbdOutNodeName << std::endl;
0453         // return Fun4AllReturnCodes::EVENT_OK;
0454     }
0455     else
0456     {
0457         MBD_charge_sum = m_mbdout->get_q(0) + m_mbdout->get_q(1);
0458     }
0459
0460     is_min_bias = (minimumbiasinfo) ? minimumbiasinfo->isAuAuMinimumBias() : false;
0461
0462     // centrality information
0463     if (!m_CentInfo)
0464     {
0465         std::cout << "SiliconSeedAnalyzer::process_event - [WARNING] - can't find CentralityInfo node " << "CentralityInfo" << std::endl;
0466         centrality_mbd_ = -1.;
0467         // return Fun4AllReturnCodes::EVENT_OK;
0468     }
0469     else
0470     {
0471         if (m_CentInfo->has_centrality_bin(CentralityInfo::PROP::mbd_NS))
0472         {
0473             centrality_mbd_ = m_CentInfo->get_centrality_bin(CentralityInfo::PROP::mbd_NS);
0474         }
0475         else
0476         {
0477             std::cout << "[WARNING/ERROR] No centrality information found in CentralityInfo. Setting centrality_mbd_ to -2. Please check!" << std::endl;
0478             m_CentInfo->identify();
0479             centrality_mbd_ = -2.;
0480         }
0481     }
0482
0483     // mbdVertex = trackerVertexX = trackerVertexY = trackerVertexZ = trackerVertexChisq = std::numeric_limits<float>::quiet_NaN();
0484     // mbdVertex = std::numeric_limits<float>::quiet_NaN();
0485     // nTracks = n_MBDVertex = n_TRKVertex = std::numeric_limits<unsigned int>::quiet_NaN();
0486     nTracks = n_MBDVertex = n_TRKVertex = 0;
0487
0488     hasMBD = false;
0489     hasTRK = false;
0490
0491     // std::cout << "Number of vertices in GlobalVertexMap: " << globalvertexmap->size() << std::endl;
0492     for (GlobalVertexMap::ConstIter iter = globalvertexmap->begin(); iter != globalvertexmap->end(); ++iter)
0493     {
0494         GlobalVertex *gvertex = iter->second;
0495
0496         if (gvertex->count_vtxs(mbdType) != 0)
0497         {
0498             // std::cout << "Found MBD vertex in GlobalVertexMap." << std::endl;
0499             hasMBD = true;
0500
0501             auto mbditer = gvertex->find_vertexes(mbdType);
0502             auto mbdvertexvector = mbditer->second;
0503
0504             n_MBDVertex += mbdvertexvector.size();
0505             for (auto &vertex : mbdvertexvector)
0506             {
0507                 MbdVertex *m_dst_vertex = m_dst_mbdvertexmap->find(vertex->get_id())->second;
0508
0509                 mbdVertexId.push_back(m_dst_vertex->get_id());
0510                 mbdVertex.push_back(m_dst_vertex->get_z());
0511                 mbdVertexCrossing.push_back(m_dst_vertex->get_beam_crossing());
0512
0513                 // std::cout << "MBD vertex z: " << m_dst_vertex->get_z() << std::endl;
0514             }
0515         }
0516
0517         if (gvertex->count_vtxs(trkType) != 0)
0518         {
0519             hasTRK = true;
0520
0521             auto trkiter = gvertex->find_vertexes(trkType);
0522             auto trkvertexvector = trkiter->second;
0523
0524             n_TRKVertex += trkvertexvector.size();
0525             for (auto &vertex : trkvertexvector)
0526             {
0527                 SvtxVertex *m_dst_vertex = m_dst_vertexmap->find(vertex->get_id())->second;
0528                 // std::cout << "Tracker vertex z: " << m_dst_vertex->get_z() << ", crossing: " << m_dst_vertex->get_beam_crossing() << std::endl;
0529
0530                 trackerVertexId.push_back(m_dst_vertex->get_id());
0531                 trackerVertexX.push_back(m_dst_vertex->get_x());
0532                 trackerVertexY.push_back(m_dst_vertex->get_y());
0533                 trackerVertexZ.push_back(m_dst_vertex->get_z());
0534                 trackerVertexChisq.push_back(m_dst_vertex->get_chisq());
0535                 trackerVertexNdof.push_back(m_dst_vertex->get_ndof());
0536                 trackerVertexNTracks.push_back(m_dst_vertex->size_tracks());
0537                 trackerVertexCrossing.push_back(m_dst_vertex->get_beam_crossing());
0538
0539                 // loop over tracks associated with this vertex and save their track IDs
0540                 std::vector<int> trackIDs;
0541                 trackIDs.clear();
0542                 for (auto trackiter = m_dst_vertex->begin_tracks(); trackiter != m_dst_vertex->end_tracks(); ++trackiter)
0543                 {
0544                     trackIDs.push_back(*trackiter);
0545                 }
0546
0547                 // print out for debugging
0548                 std::cout << "Tracker vertex ID " << m_dst_vertex->get_id() << " with crossing " << m_dst_vertex->get_beam_crossing() << " m_dst_vertex->size_tracks() = " << m_dst_vertex->size_tracks() << " trackIDs.size() = " << trackIDs.size() << ": [";
0549                 for (const auto &trackID : trackIDs)
0550                 {
0551                     std::cout << trackID << " ";
0552                 }
0553                 std::cout << "]" << std::endl;
0554
0555                 trackerVertexTrackIDs.push_back(trackIDs);
0556
0557                 // if (m_dst_vertex->get_beam_crossing() != 0)
0558                 //     continue;
0559                 // if (m_dst_vertex->size_tracks() > nTracks)
0560                 // {
0561                 //     trackerVertexX = m_dst_vertex->get_x();
0562                 //     trackerVertexY = m_dst_vertex->get_y();
0563                 //     trackerVertexZ = m_dst_vertex->get_z();
0564                 //     trackerVertexChisq = m_dst_vertex->get_chisq();
0565                 //     trackerVertexNdof = m_dst_vertex->get_ndof();
0566                 //     nTracks = m_dst_vertex->size_tracks();
0567                 // }
0568                 // if (nTracks == 0)
0569                 //     hasTRK = false;
0570             }
0571         }
0572     }
0573
0574     // loop over all vertices in MbdVertexMap and fill all vertices to ntuple
0575     // n_MBDVertex = m_dst_mbdvertexmap->size();
0576     // for (const auto& [key, vertex] : *m_dst_mbdvertexmap)
0577     // {
0578     //     mbdVertexId.push_back(vertex->get_id());
0579     //     mbdVertex.push_back(vertex->get_z());
0580     //     mbdVertexCrossing.push_back(vertex->get_beam_crossing());
0581     // }
0582
0583     // loop over all vertices in SvtxVertexMap and fill all vertices to ntuple
0584     // nSvtxVertices = m_dst_vertexmap->size();
0585     // for (const auto& [key, vertex] : *m_dst_vertexmap)
0586     // {
0587     //     trackerVertexId.push_back(vertex->get_id());
0588     //     trackerVertexX.push_back(vertex->get_x());
0589     //     trackerVertexY.push_back(vertex->get_y());
0590     //     trackerVertexZ.push_back(vertex->get_z());
0591     //     trackerVertexChisq.push_back(vertex->get_chisq()); //! For vertex from PHSimpleVertexFinder, the chisq and ndof are alway 0, need to check if this is intended
0592     //     trackerVertexNdof.push_back(vertex->get_ndof());
0593     //     trackerVertexNTracks.push_back(vertex->size_tracks());
0594     //     trackerVertexCrossing.push_back(vertex->get_beam_crossing());
0595     // }
0596
0597     // simulation setup
0598     if (isSimulation)
0599     {
0600         m_truth_info = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
0601         if (!m_truth_info)
0602         {
0603             std::cout << "[WARNING/ERROR] VertexCompare::process_event - [ERROR] - can't find G4TruthInfoContainer node " << "G4TruthInfo" << std::endl;
0604             // return Fun4AllReturnCodes::ABORTEVENT;
0605         }
0606
0607         if (!svtx_evalstack)
0608         {
0609             svtx_evalstack = new SvtxEvalStack(topNode);
0610             clustereval = svtx_evalstack->get_cluster_eval();
0611             hiteval = svtx_evalstack->get_hit_eval();
0612             truth_eval = svtx_evalstack->get_truth_eval();
0613         }
0614
0615         svtx_evalstack->next_event(topNode);
0616     }
0617
0618     FillSiliconSeedTree();
0619     FillClusterTree();
0620     if (isSimulation)
0621     {
0622         FillTruthParticleTree();
0623     }
0624
0625     outTree->Fill();
0626
0627     ++counter;
0628
0629     Cleanup();
0630
0631     return Fun4AllReturnCodes::EVENT_OK;
0632 }
0633
0634 //____________________________________________________________________________..
0635 void VertexCompare::FillSiliconSeedTree()
0636 {
0637     constexpr int kUnassociatedVertexId = -std::numeric_limits<int>::max();
0638
0639     std::vector<float> hit_x;
0640     std::vector<float> hit_y;
0641     std::vector<float> hit_z;
0642     std::vector<int> hit_rows;
0643     std::vector<int> hit_cols;
0644
0645     std::vector<int> ancestor_trackIDs;
0646     std::vector<int> ancestor_PIDs;
0647
0648     // helper function to find the vertex index based on the seed crossing and index
0649     // look up trackerVertexCrossing and trackerVertexTrackIDs, the seed crossing should match the trackerVertexCrossing
0650     // then from trackerVertexTrackIDs see which of the vertex is the seed associated with, return the vertex index
0651     // If none of the vertex's trackID list contains the seed track ID, it means the seed is not associated with any vertex
0652     // In that case, return -1 and the seed eta/phi w.r.t vertex will be set to some large or small value to indicate they are not associated with any vertex
0653     const auto find_vertex_index_for_crossing = [this](int seed_id, int crossing) -> int
0654     {
0655         // first get all vertex indices with the same crossing
0656         std::vector<int> candidate_vertex_indices;
0657         for (size_t ivtx = 0; ivtx < trackerVertexCrossing.size(); ++ivtx)
0658         {
0659             if (trackerVertexCrossing[ivtx] == crossing)
0660             {
0661                 candidate_vertex_indices.push_back(static_cast<int>(ivtx));
0662             }
0663         }
0664
0665         // then check if the seed track ID is in any of the candidate vertex track ID list
0666         for (const auto &vtx_index : candidate_vertex_indices)
0667         {
0668             const auto &trackIDs = trackerVertexTrackIDs[vtx_index];
0669             if (std::find(trackIDs.begin(), trackIDs.end(), seed_id) != trackIDs.end())
0670             {
0671                 return vtx_index;
0672             }
0673         }
0674         return -1;
0675     };
0676
0677     std::vector<std::pair<TrackSeed *, int>> seeds;
0678     std::map<int, int> first_associated_vertex_index_by_crossing;
0679     for (auto iter = silseedmap->begin(); iter != silseedmap->end(); ++iter)
0680     {
0681         auto *seed = *iter;
0682         if (!seed)
0683         {
0684             continue;
0685         }
0686
0687         const int seed_id = silseedmap->find(seed);
0688         seeds.emplace_back(seed, seed_id);
0689
0690         const int crossing = seed->get_crossing();
0691         if (first_associated_vertex_index_by_crossing.find(crossing) != first_associated_vertex_index_by_crossing.end())
0692         {
0693             continue;
0694         }
0695
0696         const int vtx_index = find_vertex_index_for_crossing(seed_id, crossing);
0697         if (vtx_index >= 0)
0698         {
0699             first_associated_vertex_index_by_crossing.emplace(crossing, vtx_index);
0700         }
0701     }
0702
0703     // loop over all silicon seeds
0704     // make a map of vector of pair of (seed id, associated vertex id) for each unique crossing for debugging purpose
0705     std::map<int, std::vector<std::pair<int, int>>> crossing_SeedIdVertexId_map;
0706     for (const auto &[seed, seed_id] : seeds)
0707     {
0708         if (!seed)
0709         {
0710             // std::cout << "VertexCompare::FillSiliconSeedTree - [ERROR] - silicon seed pointer is null, skip" << std::endl;
0711             continue;
0712         }
0713         silseed_id.push_back(seed_id);
0714         const auto si_pos = TrackSeedHelper::get_xyz(seed);
0715         silseed_x.push_back(si_pos.x());
0716         silseed_y.push_back(si_pos.y());
0717         silseed_z.push_back(si_pos.z());
0718         silseed_crossing.push_back(seed->get_crossing());
0719         if (seed->get_crossing() != SHRT_MAX)
0720         {
0721             ++nSilSeedsValidCrossing;
0722         }
0723
0724         silseed_pt.push_back(seed->get_pt());
0725         silseed_eta.push_back(seed->get_eta());
0726         silseed_phi.push_back(seed->get_phi());
0727         {
0728             const int associated_vtx_index = find_vertex_index_for_crossing(seed_id, seed->get_crossing());
0729             int eta_phi_vtx_index = associated_vtx_index;
0730             if (eta_phi_vtx_index < 0)
0731             {
0732                 const auto fallback_it = first_associated_vertex_index_by_crossing.find(seed->get_crossing());
0733                 if (fallback_it != first_associated_vertex_index_by_crossing.end())
0734                 {
0735                     eta_phi_vtx_index = fallback_it->second;
0736                 }
0737             }
0738
0739             if (associated_vtx_index >= 0)
0740             {
0741                 silseed_assocVtxId.push_back(trackerVertexId[associated_vtx_index]);
0742             }
0743             else
0744             {
0745                 silseed_assocVtxId.push_back(kUnassociatedVertexId);
0746             }
0747
0748             if (eta_phi_vtx_index >= 0)
0749             {
0750                 TVector3 seed_from_vtx(si_pos.x() - trackerVertexX[eta_phi_vtx_index], si_pos.y() - trackerVertexY[eta_phi_vtx_index], si_pos.z() - trackerVertexZ[eta_phi_vtx_index]);
0751                 silseed_eta_vtx.push_back(seed_from_vtx.Eta());
0752                 silseed_phi_vtx.push_back(seed_from_vtx.Phi());
0753             }
0754             else
0755             {
0756                 // set to some minimum value to indicate they are not associated with any vertex
0757                 silseed_eta_vtx.push_back(-1*std::numeric_limits<float>::max());
0758                 silseed_phi_vtx.push_back(-1*std::numeric_limits<float>::max());
0759             }
0760         }
0761         silseed_charge.push_back((seed->get_qOverR() > 0) ? 1 : -1);
0762         // filling crossing seed id map for debugging
0763         crossing_SeedIdVertexId_map[seed->get_crossing()].push_back(std::make_pair(seed_id, silseed_assocVtxId.back()));
0764
0765         std::vector<uint64_t> seed_cluskeys(seed->begin_cluster_keys(), seed->end_cluster_keys());
0766         silseed_clusterKeys.push_back(seed_cluskeys);
0767         silseed_cluster_layer.push_back(std::vector<unsigned int>());
0768         silseed_cluster_globalX.push_back(std::vector<float>());
0769         silseed_cluster_globalY.push_back(std::vector<float>());
0770         silseed_cluster_globalZ.push_back(std::vector<float>());
0771         silseed_cluster_phi.push_back(std::vector<float>());
0772         silseed_cluster_eta.push_back(std::vector<float>());
0773         silseed_cluster_r.push_back(std::vector<float>());
0774         silseed_cluster_phiSize.push_back(std::vector<int>());
0775         silseed_cluster_zSize.push_back(std::vector<int>());
0776         silseed_cluster_strobeID.push_back(std::vector<int>());
0777         silseed_cluster_timeBucketID.push_back(std::vector<int>());
0778         if (isSimulation)
0779         {
0780             silseed_cluster_gcluster_key.push_back(std::vector<uint64_t>());
0781             silseed_cluster_gcluster_layer.push_back(std::vector<unsigned int>());
0782             silseed_cluster_gcluster_X.push_back(std::vector<float>());
0783             silseed_cluster_gcluster_Y.push_back(std::vector<float>());
0784             silseed_cluster_gcluster_Z.push_back(std::vector<float>());
0785             silseed_cluster_gcluster_r.push_back(std::vector<float>());
0786             silseed_cluster_gcluster_phi.push_back(std::vector<float>());
0787             silseed_cluster_gcluster_eta.push_back(std::vector<float>());
0788             silseed_cluster_gcluster_edep.push_back(std::vector<float>());
0789             silseed_cluster_gcluster_adc.push_back(std::vector<int>());
0790             silseed_cluster_gcluster_phiSize.push_back(std::vector<float>());
0791             silseed_cluster_gcluster_zSize.push_back(std::vector<float>());
0792         }
0793         int nMvtx = 0, nIntt = 0;
0794         int ngMvtx = 0, ngIntt = 0;
0795         for (auto cluskey : seed_cluskeys)
0796         {
0797             auto *cluster = clustermap->findCluster(cluskey);
0798
0799             if (!cluster)
0800             {
0801                 // std::cout << "SiliconSeedAnalyzer::process_event - [ERROR] - cluster pointer is null, skip" << std::endl;
0802                 continue;
0803             }
0804
0805             unsigned int layer = TrkrDefs::getLayer(cluskey);
0806             silseed_cluster_layer[silseed_cluster_layer.size() - 1].push_back(layer);
0807
0808             auto globalpos = geometry->getGlobalPosition(cluskey, cluster);
0809             silseed_cluster_globalX[silseed_cluster_globalX.size() - 1].push_back(globalpos.x());
0810             silseed_cluster_globalY[silseed_cluster_globalY.size() - 1].push_back(globalpos.y());
0811             silseed_cluster_globalZ[silseed_cluster_globalZ.size() - 1].push_back(globalpos.z());
0812             float phi = std::atan2(globalpos.y(), globalpos.x());
0813             silseed_cluster_phi[silseed_cluster_phi.size() - 1].push_back(phi);
0814             float r = (globalpos.y() > 0) ? std::sqrt(globalpos.x() * globalpos.x() + globalpos.y() * globalpos.y()) : -std::sqrt(globalpos.x() * globalpos.x() + globalpos.y() * globalpos.y());
0815             silseed_cluster_r[silseed_cluster_r.size() - 1].push_back(r);
0816             TVector3 posvec(globalpos.x(), globalpos.y(), globalpos.z());
0817             silseed_cluster_eta[silseed_cluster_eta.size() - 1].push_back(posvec.Eta());
0818             float phisize = cluster->getPhiSize();
0819             if (phisize <= 0)
0820             {
0821                 phisize += 256;
0822             }
0823             silseed_cluster_phiSize[silseed_cluster_phiSize.size() - 1].push_back(phisize);
0824             float zsize = cluster->getZSize();
0825             silseed_cluster_zSize[silseed_cluster_zSize.size() - 1].push_back(zsize);
0826             if (TrkrDefs::getTrkrId(cluskey) == TrkrDefs::inttId)
0827             {
0828                 ++nIntt;
0829                 silseed_cluster_strobeID[silseed_cluster_strobeID.size() - 1].push_back(std::numeric_limits<int>::min());
0830                 silseed_cluster_timeBucketID[silseed_cluster_timeBucketID.size() - 1].push_back(InttDefs::getTimeBucketId(cluskey));
0831             }
0832             else if (TrkrDefs::getTrkrId(cluskey) == TrkrDefs::mvtxId)
0833             {
0834                 ++nMvtx;
0835                 silseed_cluster_strobeID[silseed_cluster_strobeID.size() - 1].push_back(MvtxDefs::getStrobeId(cluskey));
0836                 silseed_cluster_timeBucketID[silseed_cluster_timeBucketID.size() - 1].push_back(std::numeric_limits<int>::min());
0837
0838                 // save mvtx seed cluster info
0839                 mvtx_seedcluster_key.push_back(cluskey);
0840                 mvtx_seedcluster_layer.push_back(layer);
0841                 mvtx_seedcluster_chip.push_back(MvtxDefs::getChipId(cluskey));
0842                 mvtx_seedcluster_stave.push_back(MvtxDefs::getStaveId(cluskey));
0843                 mvtx_seedcluster_globalX.push_back(globalpos.x());
0844                 mvtx_seedcluster_globalY.push_back(globalpos.y());
0845                 mvtx_seedcluster_globalZ.push_back(globalpos.z());
0846                 mvtx_seedcluster_phi.push_back(phi);
0847                 mvtx_seedcluster_eta.push_back(posvec.Eta());
0848                 mvtx_seedcluster_r.push_back(r);
0849                 mvtx_seedcluster_phiSize.push_back(phisize);
0850                 mvtx_seedcluster_zSize.push_back(zsize);
0851                 mvtx_seedcluster_strobeID.push_back(MvtxDefs::getStrobeId(cluskey));
0852                 mvtx_seedcluster_matchedcrossing.push_back(seed->get_crossing());
0853                 // get hitrow and hitcol from cluster hit assoc
0854                 {
0855                     Clean(hit_x);
0856                     Clean(hit_y);
0857                     Clean(hit_z);
0858                     Clean(hit_rows);
0859                     Clean(hit_cols);
0860
0861                     TrkrClusterHitAssoc::ConstRange hitrange = clusterhitassoc->getHits(cluskey);
0862                     for (TrkrClusterHitAssoc::ConstIterator hititer = hitrange.first; hititer != hitrange.second; ++hititer)
0863                     {
0864                         TrkrDefs::hitkey hitkey = hititer->second;
0865
0866                         int hitrow = MvtxDefs::getRow(hitkey);
0867                         int hitcol = MvtxDefs::getCol(hitkey);
0868                         hit_rows.push_back(hitrow);
0869                         hit_cols.push_back(hitcol);
0870
0871                         MvtxPixelDefs::pixelkey pixelkey = MvtxPixelDefs::gen_pixelkey_from_coors(layer, MvtxDefs::getStaveId(cluskey), MvtxDefs::getChipId(cluskey), hitrow, hitcol);
0872                         uint32_t hitsetkey = MvtxPixelDefs::get_hitsetkey(pixelkey);
0873
0874                         float localX = std::numeric_limits<float>::quiet_NaN();
0875                         float localZ = std::numeric_limits<float>::quiet_NaN();
0876                         SegmentationAlpide::detectorToLocal(hitrow, hitcol, localX, localZ);
0877                         Acts::Vector2 local(localX * Acts::UnitConstants::cm, localZ * Acts::UnitConstants::cm);
0878
0879                         const auto &surface = geometry->maps().getSiliconSurface(hitsetkey);
0880                         auto glob = surface->localToGlobal(geometry->geometry().getGeoContext(), local, Acts::Vector3());
0881
0882                         hit_x.push_back(glob.x() / Acts::UnitConstants::cm);
0883                         hit_y.push_back(glob.y() / Acts::UnitConstants::cm);
0884                         hit_z.push_back(glob.z() / Acts::UnitConstants::cm);
0885                     }
0886
0887                     mvtx_seedcluster_hitX.push_back(hit_x);
0888                     mvtx_seedcluster_hitY.push_back(hit_y);
0889                     mvtx_seedcluster_hitZ.push_back(hit_z);
0890                     mvtx_seedcluster_hitrow.push_back(hit_rows);
0891                     mvtx_seedcluster_hitcol.push_back(hit_cols);
0892                 }
0893
0894                 // if simulation, get matched G4 particle info
0895                 if (isSimulation)
0896                 {
0897                     Clean(ancestor_trackIDs);
0898                     Clean(ancestor_PIDs);
0899
0900                     PHG4Particle *ptcl = clustereval->max_truth_particle_by_energy(cluskey); // alternatively, can also try max_truth_particle_by_cluster_energy
0901                     if (ptcl)
0902                     {
0903                         mvtx_seedcluster_matchedG4P_trackID.push_back(ptcl->get_track_id());
0904                         mvtx_seedcluster_matchedG4P_PID.push_back(ptcl->get_pid());
0905                         mvtx_seedcluster_matchedG4P_E.push_back(ptcl->get_e());
0906                         ROOT::Math::PtEtaPhiEVector g4p4(ptcl->get_px(), ptcl->get_py(), ptcl->get_pz(), ptcl->get_e());
0907                         mvtx_seedcluster_matchedG4P_pT.push_back(g4p4.Pt());
0908                         mvtx_seedcluster_matchedG4P_eta.push_back(g4p4.Eta());
0909                         mvtx_seedcluster_matchedG4P_phi.push_back(g4p4.Phi());
0910
0911                         // get ancestor info
0912                         PHG4Particle *ancestor = m_truth_info->GetParticle(ptcl->get_parent_id());
0913                         while (ancestor != nullptr)
0914                         {
0915                             ancestor_trackIDs.push_back(ancestor->get_track_id());
0916                             ancestor_PIDs.push_back(ancestor->get_pid());
0917                             ancestor = m_truth_info->GetParticle(ancestor->get_parent_id());
0918                         }
0919                         mvtx_seedcluster_matchedG4P_ancestor_trackID.push_back(ancestor_trackIDs);
0920                         mvtx_seedcluster_matchedG4P_ancestor_PID.push_back(ancestor_PIDs);
0921                     }
0922                     else
0923                     {
0924                         std::cout << "VertexCompare::FillSiliconSeedTree - [WARNING] - no matched G4 particle found for cluster " << cluskey << std::endl;
0925                         mvtx_seedcluster_matchedG4P_trackID.push_back(std::numeric_limits<int>::max());
0926                         mvtx_seedcluster_matchedG4P_PID.push_back(std::numeric_limits<int>::max());
0927                         mvtx_seedcluster_matchedG4P_E.push_back(std::numeric_limits<float>::quiet_NaN());
0928                         mvtx_seedcluster_matchedG4P_pT.push_back(std::numeric_limits<float>::quiet_NaN());
0929                         mvtx_seedcluster_matchedG4P_eta.push_back(std::numeric_limits<float>::quiet_NaN());
0930                         mvtx_seedcluster_matchedG4P_phi.push_back(std::numeric_limits<float>::quiet_NaN());
0931                         mvtx_seedcluster_matchedG4P_ancestor_trackID.push_back(std::vector<int>());
0932                         mvtx_seedcluster_matchedG4P_ancestor_PID.push_back(std::vector<int>());
0933                     }
0934
0935                 }
0936             }
0937             else
0938             {
0939                 std::cout << "SiliconSeedAnalyzer::process_event - [WARNING] - cluster is neither INTT nor MVTX. Please check." << std::endl;
0940                 silseed_cluster_strobeID[silseed_cluster_strobeID.size() - 1].push_back(std::numeric_limits<int>::max());
0941                 silseed_cluster_timeBucketID[silseed_cluster_timeBucketID.size() - 1].push_back(std::numeric_limits<int>::max());
0942             }
0943
0944             if (isSimulation)
0945             {
0946                 std::pair<TrkrDefs::cluskey, std::shared_ptr<TrkrCluster>> truthclus = clustereval->max_truth_cluster_by_energy(cluskey);
0947                 const auto truth_key = truthclus.first;
0948                 const auto &truth_cluster = truthclus.second;
0949                 if (truth_cluster)
0950                 {
0951                     const float gx = truth_cluster->getX();
0952                     const float gy = truth_cluster->getY();
0953                     const float gz = truth_cluster->getZ();
0954                     TVector3 gpos(gx, gy, gz);
0955
0956                     silseed_cluster_gcluster_key.back().push_back(truth_key);
0957                     silseed_cluster_gcluster_layer.back().push_back(TrkrDefs::getLayer(truth_key));
0958                     silseed_cluster_gcluster_X.back().push_back(gx);
0959                     silseed_cluster_gcluster_Y.back().push_back(gy);
0960                     silseed_cluster_gcluster_Z.back().push_back(gz);
0961                     silseed_cluster_gcluster_r.back().push_back(std::sqrt(gx * gx + gy * gy));
0962                     silseed_cluster_gcluster_phi.back().push_back(gpos.Phi());
0963                     silseed_cluster_gcluster_eta.back().push_back(gpos.Eta());
0964                     silseed_cluster_gcluster_edep.back().push_back(truth_cluster->getError(0, 0));
0965                     silseed_cluster_gcluster_adc.back().push_back(truth_cluster->getAdc());
0966                     silseed_cluster_gcluster_phiSize.back().push_back(truth_cluster->getSize(1, 1));
0967                     silseed_cluster_gcluster_zSize.back().push_back(truth_cluster->getSize(2, 2));
0968
0969                     if (TrkrDefs::getTrkrId(cluskey) == TrkrDefs::mvtxId)
0970                     {
0971                         ++ngMvtx;
0972                     }
0973                     else if (TrkrDefs::getTrkrId(cluskey) == TrkrDefs::inttId)
0974                     {
0975                         ++ngIntt;
0976                     }
0977                 }
0978                 else
0979                 {
0980                     silseed_cluster_gcluster_key.back().push_back(0);
0981                     silseed_cluster_gcluster_layer.back().push_back(std::numeric_limits<unsigned int>::max());
0982                     silseed_cluster_gcluster_X.back().push_back(-1 * std::numeric_limits<float>::max());
0983                     silseed_cluster_gcluster_Y.back().push_back(-1 * std::numeric_limits<float>::max());
0984                     silseed_cluster_gcluster_Z.back().push_back(-1 * std::numeric_limits<float>::max());
0985                     silseed_cluster_gcluster_r.back().push_back(-1 * std::numeric_limits<float>::max());
0986                     silseed_cluster_gcluster_phi.back().push_back(-1 * std::numeric_limits<float>::max());
0987                     silseed_cluster_gcluster_eta.back().push_back(-1 * std::numeric_limits<float>::max());
0988                     silseed_cluster_gcluster_edep.back().push_back(-1 * std::numeric_limits<float>::max());
0989                     silseed_cluster_gcluster_adc.back().push_back(-1 * std::numeric_limits<int>::max());
0990                     silseed_cluster_gcluster_phiSize.back().push_back(-1 * std::numeric_limits<float>::max());
0991                     silseed_cluster_gcluster_zSize.back().push_back(-1 * std::numeric_limits<float>::max());
0992                 }
0993             }
0994         }
0995         silseed_nMvtx.push_back(nMvtx);
0996         silseed_nIntt.push_back(nIntt);
0997         if (isSimulation)
0998         {
0999             silseed_ngmvtx.push_back(ngMvtx);
1000             silseed_ngintt.push_back(ngIntt);
1001         }
1002     }
1003     nTotalSilSeeds = silseed_id.size();
1004
1005     std::cout << "Total silicon seeds in this event: " << nTotalSilSeeds << " with valid crossing: " << nSilSeedsValidCrossing << std::endl;
1006
1007
1008     // now print out the crossing seed id map for debugging
1009     std::cout << "Crossing to seed ID mapping:" << std::endl;
1010     for (const auto &[crossing, seed_id_vertex_id_pairs] : crossing_SeedIdVertexId_map)
1011     {
1012         std::cout << "  Crossing " << crossing << ": Seed IDs [";
1013         for (size_t i = 0; i < seed_id_vertex_id_pairs.size(); ++i)
1014         {
1015             std::cout << "(" << seed_id_vertex_id_pairs[i].first << ", " << seed_id_vertex_id_pairs[i].second << ")";
1016             if (i != seed_id_vertex_id_pairs.size() - 1)
1017             {
1018                 std::cout << ", ";
1019             }
1020         }
1021         std::cout << "]" << std::endl;
1022     }
1023
1024 }
1025
1026 //____________________________________________________________________________..
1027 void VertexCompare::FillClusterTree()
1028 {
1029     for (const auto &det : {TrkrDefs::TrkrId::inttId})
1030     {
1031         for (const auto &hitsetkey : clustermap->getHitSetKeys(det))
1032         {
1033             auto range = clustermap->getClusters(hitsetkey);
1034             for (auto iter = range.first; iter != range.second; ++iter)
1035             {
1036                 uint64_t key = iter->first;
1037                 clusterKey.push_back(key);
1038                 unsigned int layer = TrkrDefs::getLayer(key);
1039                 cluster_layer.push_back(layer);
1040                 auto *cluster = iter->second;
1041                 auto globalpos = geometry->getGlobalPosition(key, cluster);
1042                 cluster_globalX.push_back(globalpos.x());
1043                 cluster_globalY.push_back(globalpos.y());
1044                 cluster_globalZ.push_back(globalpos.z());
1045                 cluster_phi.push_back(std::atan2(globalpos.y(), globalpos.x()));
1046                 TVector3 posvec(globalpos.x(), globalpos.y(), globalpos.z());
1047                 cluster_eta.push_back(posvec.Eta());
1048                 cluster_r.push_back((globalpos.y() > 0) ? std::sqrt(globalpos.x() * globalpos.x() + globalpos.y() * globalpos.y()) : -std::sqrt(globalpos.x() * globalpos.x() + globalpos.y() * globalpos.y()));
1049                 int phiSize = cluster->getPhiSize();
1050                 if (phiSize <= 0)
1051                 {
1052                     phiSize += 256;
1053                 }
1054                 cluster_phiSize.push_back(phiSize);
1055                 int zSize = cluster->getZSize();
1056                 cluster_zSize.push_back(zSize);
1057                 cluster_adc.push_back(cluster->getAdc());
1058                 switch (det)
1059                 {
1060                 case TrkrDefs::TrkrId::mvtxId:
1061                     cluster_chip.push_back(MvtxDefs::getChipId(key));
1062                     cluster_stave.push_back(MvtxDefs::getStaveId(key));
1063                     cluster_timeBucketID.push_back(std::numeric_limits<int>::min());
1064                     cluster_LocalX.push_back(cluster->getLocalX());
1065                     cluster_LocalY.push_back(cluster->getLocalY());
1066                     break;
1067                 case TrkrDefs::TrkrId::inttId:
1068                     cluster_ladderZId.push_back(InttDefs::getLadderZId(key));
1069                     cluster_ladderPhiId.push_back(InttDefs::getLadderPhiId(key));
1070                     cluster_timeBucketID.push_back(InttDefs::getTimeBucketId(key));
1071                     cluster_LocalX.push_back(cluster->getLocalX());
1072                     cluster_LocalY.push_back(cluster->getLocalY());
1073                     break;
1074                 default:
1075                     std::cout << "SiliconSeedAnalyzer::process_event - [WARNING] - cluster is neither INTT nor MVTX. Please check." << std::endl;
1076                     cluster_chip.push_back(std::numeric_limits<int>::max());
1077                     cluster_stave.push_back(std::numeric_limits<int>::max());
1078                     cluster_ladderZId.push_back(std::numeric_limits<int>::max());
1079                     cluster_ladderPhiId.push_back(std::numeric_limits<int>::max());
1080                     cluster_timeBucketID.push_back(std::numeric_limits<int>::max());
1081                     cluster_LocalX.push_back(std::numeric_limits<float>::max());
1082                     cluster_LocalY.push_back(std::numeric_limits<float>::max());
1083                     break;
1084                 }
1085
1086                 if (isSimulation)
1087                 {
1088                     PHG4Particle *ptcl_maxE = (clustereval) ? clustereval->max_truth_particle_by_energy(key) : nullptr;
1089                     if (ptcl_maxE)
1090                     {
1091                         cluster_matchedG4P_trackID.push_back(ptcl_maxE->get_track_id());
1092                         cluster_matchedG4P_PID.push_back(ptcl_maxE->get_pid());
1093                         cluster_matchedG4P_E.push_back(ptcl_maxE->get_e());
1094
1095                         ROOT::Math::PxPyPzEVector g4p4(ptcl_maxE->get_px(), ptcl_maxE->get_py(), ptcl_maxE->get_pz(), ptcl_maxE->get_e());
1096                         cluster_matchedG4P_pT.push_back(g4p4.Pt());
1097                         cluster_matchedG4P_eta.push_back(g4p4.Eta());
1098                         cluster_matchedG4P_phi.push_back(g4p4.Phi());
1099                     }
1100                     else
1101                     {
1102                         cluster_matchedG4P_trackID.push_back(std::numeric_limits<int>::max());
1103                         cluster_matchedG4P_PID.push_back(std::numeric_limits<int>::max());
1104                         cluster_matchedG4P_E.push_back(-1*std::numeric_limits<float>::max());
1105                         cluster_matchedG4P_pT.push_back(-1*std::numeric_limits<float>::max());
1106                         cluster_matchedG4P_eta.push_back(-1*std::numeric_limits<float>::max());
1107                         cluster_matchedG4P_phi.push_back(-1*std::numeric_limits<float>::max());
1108                     }
1109                 }
1110             }
1111         }
1112     }
1113 }
1114
1115 //____________________________________________________________________________..
1116 void VertexCompare::FillTruthParticleTree()
1117 {
1118     if (!m_truth_info)
1119     {
1120         std::cout << "VertexCompare::FillTruthParticleTree - [WARNING] - missing G4TruthInfo, skip filling truth particle branches" << std::endl;
1121         return;
1122     }
1123
1124     // truth primary vertices
1125     const auto vtx_range = m_truth_info->GetPrimaryVtxRange();
1126     for (auto iter = vtx_range.first; iter != vtx_range.second; ++iter)
1127     {
1128         const int point_id = iter->first;
1129         PHG4VtxPoint *point = iter->second;
1130         if (!point)
1131         {
1132             continue;
1133         }
1134
1135         ++nTruthVertex;
1136         if (m_truth_info->isEmbededVtx(point_id) == 0)
1137         {
1138             TruthVertexX.push_back(point->get_x());
1139             TruthVertexY.push_back(point->get_y());
1140             TruthVertexZ.push_back(point->get_z());
1141         }
1142     }
1143
1144     auto fill_ancestor_info = [this](PHG4Particle *ptcl, std::vector<int> &ancestor_trackIDs, std::vector<int> &ancestor_PIDs)
1145     {
1146         PHG4Particle *ancestor = m_truth_info->GetParticle(ptcl->get_parent_id());
1147         while (ancestor != nullptr)
1148         {
1149             ancestor_trackIDs.push_back(ancestor->get_track_id());
1150             ancestor_PIDs.push_back(ancestor->get_pid());
1151             ancestor = m_truth_info->GetParticle(ancestor->get_parent_id());
1152         }
1153     };
1154
1155     auto fill_particle_kinematics = [](PHG4Particle *ptcl, std::vector<float> &out_pT, std::vector<float> &out_eta, std::vector<float> &out_phi, std::vector<float> &out_E, std::vector<int> &out_PID, std::vector<int> &out_trackID)
1156     {
1157         ROOT::Math::PxPyPzEVector p4(ptcl->get_px(), ptcl->get_py(), ptcl->get_pz(), ptcl->get_e());
1158         out_pT.push_back(p4.Pt());
1159         out_eta.push_back(p4.Eta());
1160         out_phi.push_back(p4.Phi());
1161         out_E.push_back(ptcl->get_e());
1162         out_PID.push_back(ptcl->get_pid());
1163         out_trackID.push_back(ptcl->get_track_id());
1164     };
1165
1166     auto fill_primary_particle_info = [](PHG4Particle *ptcl, std::vector<TString> &out_particleClass, std::vector<bool> &out_isStable, std::vector<double> &out_charge, std::vector<bool> &out_isChargedHadron)
1167     {
1168         TString particleClass = "Unknown";
1169         bool isStable = false;
1170         double charge = 0;
1171
1172         TParticlePDG *pdg = TDatabasePDG::Instance()->GetParticle(ptcl->get_pid());
1173         if (pdg)
1174         {
1175             particleClass = TString(pdg->ParticleClass());
1176             isStable = (pdg->Stable() == 1);
1177             charge = pdg->Charge();
1178         }
1179
1180         bool isHadron = (particleClass.Contains("Baryon") || particleClass.Contains("Meson"));
1181         bool isChargedHadron = (isStable && (charge != 0) && isHadron);
1182
1183         out_particleClass.push_back(particleClass);
1184         out_isStable.push_back(isStable);
1185         out_charge.push_back(charge);
1186         out_isChargedHadron.push_back(isChargedHadron);
1187     };
1188
1189     auto fill_truthcluster_match_info = [this](PHG4Particle *ptcl,                                        //
1190                                                std::vector<std::vector<float>> &out_truthcluster_X,       //
1191                                                std::vector<std::vector<float>> &out_truthcluster_Y,       //
1192                                                std::vector<std::vector<float>> &out_truthcluster_Z,       //
1193                                                std::vector<std::vector<float>> &out_truthcluster_edep,    //
1194                                                std::vector<std::vector<float>> &out_truthcluster_adc,     //
1195                                                std::vector<std::vector<float>> &out_truthcluster_r,       //
1196                                                std::vector<std::vector<float>> &out_truthcluster_phi,     //
1197                                                std::vector<std::vector<float>> &out_truthcluster_eta,     //
1198                                                std::vector<std::vector<float>> &out_truthcluster_phisize, //
1199                                                std::vector<std::vector<float>> &out_truthcluster_zsize,   //
1200                                                std::vector<std::vector<float>> &out_recocluster_globalX,  //
1201                                                std::vector<std::vector<float>> &out_recocluster_globalY,  //
1202                                                std::vector<std::vector<float>> &out_recocluster_globalZ,  //
1203                                                std::vector<std::vector<float>> &out_recocluster_r,        //
1204                                                std::vector<std::vector<float>> &out_recocluster_phi,      //
1205                                                std::vector<std::vector<float>> &out_recocluster_eta,      //
1206                                                std::vector<std::vector<float>> &out_recocluster_phisize,  //
1207                                                std::vector<std::vector<float>> &out_recocluster_zsize,    //
1208                                                std::vector<std::vector<float>> &out_recocluster_adc       //
1209                                         )
1210     {
1211         std::vector<float> truthcluster_X;
1212         std::vector<float> truthcluster_Y;
1213         std::vector<float> truthcluster_Z;
1214         std::vector<float> truthcluster_edep;
1215         std::vector<float> truthcluster_adc;
1216         std::vector<float> truthcluster_r;
1217         std::vector<float> truthcluster_phi;
1218         std::vector<float> truthcluster_eta;
1219         std::vector<float> truthcluster_phisize;
1220         std::vector<float> truthcluster_zsize;
1221         std::vector<float> recocluster_globalX;
1222         std::vector<float> recocluster_globalY;
1223         std::vector<float> recocluster_globalZ;
1224         std::vector<float> recocluster_r;
1225         std::vector<float> recocluster_phi;
1226         std::vector<float> recocluster_eta;
1227         std::vector<float> recocluster_phisize;
1228         std::vector<float> recocluster_zsize;
1229         std::vector<float> recocluster_adc;
1230
1231         if (truth_eval && clustereval)
1232         {
1233             const auto truth_clusters = truth_eval->all_truth_clusters(ptcl);
1234             for (const auto &[ckey, gclus] : truth_clusters)
1235             {
1236                 if (!gclus)
1237                 {
1238                     continue;
1239                 }
1240
1241                 // only get cluster in MVTX and INTT for now and skip TPC
1242                 if (TrkrDefs::getTrkrId(ckey) != TrkrDefs::mvtxId && TrkrDefs::getTrkrId(ckey) != TrkrDefs::inttId)
1243                 {
1244                     continue;
1245                 }
1246
1247                 const float gx = gclus->getX();
1248                 const float gy = gclus->getY();
1249                 const float gz = gclus->getZ();
1250                 const float gedep = gclus->getError(0, 0);
1251                 const float gadc = static_cast<float>(gclus->getAdc());
1252
1253                 TVector3 gpos(gx, gy, gz);
1254                 truthcluster_X.push_back(gx);
1255                 truthcluster_Y.push_back(gy);
1256                 truthcluster_Z.push_back(gz);
1257                 truthcluster_edep.push_back(gedep);
1258                 truthcluster_adc.push_back(gadc);
1259                 truthcluster_r.push_back(std::sqrt(gx * gx + gy * gy));
1260                 truthcluster_phi.push_back(gpos.Phi());
1261                 truthcluster_eta.push_back(gpos.Eta());
1262                 truthcluster_phisize.push_back(gclus->getSize(1, 1));
1263                 truthcluster_zsize.push_back(gclus->getSize(2, 2));
1264
1265                 float x = std::numeric_limits<float>::quiet_NaN();
1266                 float y = std::numeric_limits<float>::quiet_NaN();
1267                 float z = std::numeric_limits<float>::quiet_NaN();
1268                 float r = std::numeric_limits<float>::quiet_NaN();
1269                 float phi = std::numeric_limits<float>::quiet_NaN();
1270                 float eta = std::numeric_limits<float>::quiet_NaN();
1271                 float phisize = std::numeric_limits<float>::quiet_NaN();
1272                 float zsize = std::numeric_limits<float>::quiet_NaN();
1273                 float adc = std::numeric_limits<float>::quiet_NaN();
1274
1275                 const auto [reco_ckey, reco_cluster] = clustereval->reco_cluster_from_truth_cluster(ckey, gclus);
1276                 if (reco_cluster)
1277                 {
1278                     const auto global = geometry->getGlobalPosition(reco_ckey, reco_cluster);
1279                     x = global.x();
1280                     y = global.y();
1281                     z = global.z();
1282
1283                     TVector3 reco_pos(x, y, z);
1284                     r = std::sqrt(x * x + y * y);
1285                     phi = reco_pos.Phi();
1286                     eta = reco_pos.Eta();
1287                     phisize = reco_cluster->getPhiSize();
1288                     zsize = reco_cluster->getZSize();
1289                     adc = static_cast<float>(reco_cluster->getAdc());
1290                 }
1291
1292                 recocluster_globalX.push_back(x);
1293                 recocluster_globalY.push_back(y);
1294                 recocluster_globalZ.push_back(z);
1295                 recocluster_r.push_back(r);
1296                 recocluster_phi.push_back(phi);
1297                 recocluster_eta.push_back(eta);
1298                 recocluster_phisize.push_back(phisize);
1299                 recocluster_zsize.push_back(zsize);
1300                 recocluster_adc.push_back(adc);
1301             }
1302         }
1303
1304         out_truthcluster_X.push_back(truthcluster_X);
1305         out_truthcluster_Y.push_back(truthcluster_Y);
1306         out_truthcluster_Z.push_back(truthcluster_Z);
1307         out_truthcluster_edep.push_back(truthcluster_edep);
1308         out_truthcluster_adc.push_back(truthcluster_adc);
1309         out_truthcluster_r.push_back(truthcluster_r);
1310         out_truthcluster_phi.push_back(truthcluster_phi);
1311         out_truthcluster_eta.push_back(truthcluster_eta);
1312         out_truthcluster_phisize.push_back(truthcluster_phisize);
1313         out_truthcluster_zsize.push_back(truthcluster_zsize);
1314         out_recocluster_globalX.push_back(recocluster_globalX);
1315         out_recocluster_globalY.push_back(recocluster_globalY);
1316         out_recocluster_globalZ.push_back(recocluster_globalZ);
1317         out_recocluster_r.push_back(recocluster_r);
1318         out_recocluster_phi.push_back(recocluster_phi);
1319         out_recocluster_eta.push_back(recocluster_eta);
1320         out_recocluster_phisize.push_back(recocluster_phisize);
1321         out_recocluster_zsize.push_back(recocluster_zsize);
1322         out_recocluster_adc.push_back(recocluster_adc);
1323     };
1324
1325     std::vector<int> ancestor_trackIDs;
1326     std::vector<int> ancestor_PIDs;
1327
1328     // primary PHG4 particles
1329     const auto primary_range = m_truth_info->GetPrimaryParticleRange();
1330     for (auto iter = primary_range.first; iter != primary_range.second; ++iter)
1331     {
1332         PHG4Particle *ptcl = iter->second;
1333         if (!ptcl)
1334         {
1335             continue;
1336         }
1337
1338         Clean(ancestor_trackIDs);
1339         Clean(ancestor_PIDs);
1340         fill_particle_kinematics(ptcl, PrimaryPHG4Ptcl_pT, PrimaryPHG4Ptcl_eta, PrimaryPHG4Ptcl_phi, PrimaryPHG4Ptcl_E, PrimaryPHG4Ptcl_PID, PrimaryPHG4Ptcl_trackID);
1341         fill_primary_particle_info(ptcl, PrimaryPHG4Ptcl_ParticleClass, PrimaryPHG4Ptcl_isStable, PrimaryPHG4Ptcl_charge, PrimaryPHG4Ptcl_isChargedHadron);
1342         fill_ancestor_info(ptcl, ancestor_trackIDs, ancestor_PIDs);
1343         PrimaryPHG4Ptcl_ancestor_trackID.push_back(ancestor_trackIDs);
1344         PrimaryPHG4Ptcl_ancestor_PID.push_back(ancestor_PIDs);
1345         fill_truthcluster_match_info(ptcl, PrimaryPHG4Ptcl_truthcluster_X, PrimaryPHG4Ptcl_truthcluster_Y, PrimaryPHG4Ptcl_truthcluster_Z, PrimaryPHG4Ptcl_truthcluster_edep, PrimaryPHG4Ptcl_truthcluster_adc, PrimaryPHG4Ptcl_truthcluster_r, PrimaryPHG4Ptcl_truthcluster_phi, PrimaryPHG4Ptcl_truthcluster_eta, PrimaryPHG4Ptcl_truthcluster_phisize, PrimaryPHG4Ptcl_truthcluster_zsize,
1346                                      PrimaryPHG4Ptcl_recocluster_globalX, PrimaryPHG4Ptcl_recocluster_globalY, PrimaryPHG4Ptcl_recocluster_globalZ, PrimaryPHG4Ptcl_recocluster_r, PrimaryPHG4Ptcl_recocluster_phi, PrimaryPHG4Ptcl_recocluster_eta, PrimaryPHG4Ptcl_recocluster_phisize, PrimaryPHG4Ptcl_recocluster_zsize, PrimaryPHG4Ptcl_recocluster_adc);
1347     }
1348     // clean ancestor info vectors before reusing for sPHENIX primary particles
1349     Clean(ancestor_trackIDs);
1350     Clean(ancestor_PIDs);
1351
1352     const auto sPHENIXprimary_particle_range = m_truth_info->GetSPHENIXPrimaryParticleRange();
1353     std::cout << "Number of sPHENIX primary particles: " << std::distance(sPHENIXprimary_particle_range.first, sPHENIXprimary_particle_range.second) << std::endl;
1354     for (auto iter = sPHENIXprimary_particle_range.first; iter != sPHENIXprimary_particle_range.second; ++iter)
1355     {
1356         PHG4Particle *ptcl = iter->second;
1357         if (!ptcl)
1358         {
1359             continue;
1360         }
1361         // get the particle in truth info container with the same track ID
1362         PHG4Particle *real_ptcl = m_truth_info->GetParticle(ptcl->get_track_id());
1363
1364         Clean(ancestor_trackIDs);
1365         Clean(ancestor_PIDs);
1366         fill_particle_kinematics(real_ptcl, sPHENIXPrimary_pT, sPHENIXPrimary_eta, sPHENIXPrimary_phi, sPHENIXPrimary_E, sPHENIXPrimary_PID, sPHENIXPrimary_trackID);
1367         fill_primary_particle_info(real_ptcl, sPHENIXPrimary_ParticleClass, sPHENIXPrimary_isStable, sPHENIXPrimary_charge, sPHENIXPrimary_isChargedHadron);
1368         fill_ancestor_info(real_ptcl, ancestor_trackIDs, ancestor_PIDs);
1369         sPHENIXPrimary_ancestor_trackID.push_back(ancestor_trackIDs);
1370         sPHENIXPrimary_ancestor_PID.push_back(ancestor_PIDs);
1371         fill_truthcluster_match_info(real_ptcl, sPHENIXPrimary_truthcluster_X, sPHENIXPrimary_truthcluster_Y, sPHENIXPrimary_truthcluster_Z, sPHENIXPrimary_truthcluster_edep, sPHENIXPrimary_truthcluster_adc, sPHENIXPrimary_truthcluster_r, sPHENIXPrimary_truthcluster_phi, sPHENIXPrimary_truthcluster_eta, sPHENIXPrimary_truthcluster_phisize, sPHENIXPrimary_truthcluster_zsize,
1372                                      sPHENIXPrimary_recocluster_globalX, sPHENIXPrimary_recocluster_globalY, sPHENIXPrimary_recocluster_globalZ, sPHENIXPrimary_recocluster_r, sPHENIXPrimary_recocluster_phi, sPHENIXPrimary_recocluster_eta, sPHENIXPrimary_recocluster_phisize, sPHENIXPrimary_recocluster_zsize, sPHENIXPrimary_recocluster_adc);
1373
1374         // print out the truth particle info and how many truth and reco clusters are matched for debugging
1375         std::cout << "sPHENIX Primary Particle - trackID: " << real_ptcl->get_track_id() << ", PID: " << real_ptcl->get_pid() << ", pT: " << sPHENIXPrimary_pT.back() << ", eta: " << sPHENIXPrimary_eta.back() << ", phi: " << sPHENIXPrimary_phi.back() << ", E: " << sPHENIXPrimary_E.back() << std::endl;
1376         std::cout << "  Matched truth clusters: " << sPHENIXPrimary_truthcluster_X.back().size() << std::endl;
1377         std::cout << "  Matched reco clusters: " << sPHENIXPrimary_recocluster_globalX.back().size() << std::endl;
1378         // flag a particle if it has more reco clusters than truth clusters
1379         if (sPHENIXPrimary_recocluster_globalX.back().size() > sPHENIXPrimary_truthcluster_X.back().size())
1380         {
1381             std::cout << "  *** Particle has more reco clusters than truth clusters ***" << std::endl;
1382         }
1383     }
1384     // again clean it before reusing for all PHG4 particles
1385     Clean(ancestor_trackIDs);
1386     Clean(ancestor_PIDs);
1387
1388     // all PHG4 particles
1389     const auto all_particle_range = m_truth_info->GetParticleRange();
1390     for (auto iter = all_particle_range.first; iter != all_particle_range.second; ++iter)
1391     {
1392         PHG4Particle *ptcl = iter->second;
1393         if (!ptcl)
1394         {
1395             continue;
1396         }
1397
1398         Clean(ancestor_trackIDs);
1399         Clean(ancestor_PIDs);
1400         fill_particle_kinematics(ptcl, AllPHG4Ptcl_pT, AllPHG4Ptcl_eta, AllPHG4Ptcl_phi, AllPHG4Ptcl_E, AllPHG4Ptcl_PID, AllPHG4Ptcl_trackID);
1401         fill_ancestor_info(ptcl, ancestor_trackIDs, ancestor_PIDs);
1402         AllPHG4Ptcl_ancestor_trackID.push_back(ancestor_trackIDs);
1403         AllPHG4Ptcl_ancestor_PID.push_back(ancestor_PIDs);
1404         fill_truthcluster_match_info(ptcl, AllPHG4Ptcl_truthcluster_X, AllPHG4Ptcl_truthcluster_Y, AllPHG4Ptcl_truthcluster_Z, AllPHG4Ptcl_truthcluster_edep, AllPHG4Ptcl_truthcluster_adc, AllPHG4Ptcl_truthcluster_r, AllPHG4Ptcl_truthcluster_phi, AllPHG4Ptcl_truthcluster_eta, AllPHG4Ptcl_truthcluster_phisize, AllPHG4Ptcl_truthcluster_zsize, AllPHG4Ptcl_recocluster_globalX,
1405                                      AllPHG4Ptcl_recocluster_globalY, AllPHG4Ptcl_recocluster_globalZ, AllPHG4Ptcl_recocluster_r, AllPHG4Ptcl_recocluster_phi, AllPHG4Ptcl_recocluster_eta, AllPHG4Ptcl_recocluster_phisize, AllPHG4Ptcl_recocluster_zsize, AllPHG4Ptcl_recocluster_adc);
1406     }
1407
1408     N_PrimaryPHG4Ptcl = PrimaryPHG4Ptcl_trackID.size();
1409     N_sPHENIXPrimary = sPHENIXPrimary_trackID.size();
1410     N_AllPHG4Ptcl = AllPHG4Ptcl_trackID.size();
1411     // final clean up of ancestor info vectors
1412     Clean(ancestor_trackIDs);
1413     Clean(ancestor_PIDs);
1414 }
1415
1416 //____________________________________________________________________________..
1417 void VertexCompare::Cleanup()
1418 {
1419     nTotalSilSeeds = 0;
1420     nSilSeedsValidCrossing = 0;
1421     MBD_charge_sum = 0;
1422     N_PrimaryPHG4Ptcl = 0;
1423     N_sPHENIXPrimary = 0;
1424     N_AllPHG4Ptcl = 0;
1425     nTruthVertex = 0;
1426
1427     Clean(firedTriggers);
1428     Clean(TruthVertexX);
1429     Clean(TruthVertexY);
1430     Clean(TruthVertexZ);
1431
1432     Clean(mbdVertex);
1433     Clean(mbdVertexId);
1434     Clean(mbdVertexCrossing);
1435
1436     Clean(trackerVertexId);
1437     Clean(trackerVertexX);
1438     Clean(trackerVertexY);
1439     Clean(trackerVertexZ);
1440     Clean(trackerVertexChisq);
1441     Clean(trackerVertexNdof);
1442     Clean(trackerVertexNTracks);
1443     Clean(trackerVertexCrossing);
1444     Clean(trackerVertexTrackIDs);
1445
1446     Clean(silseed_id);
1447     Clean(silseed_assocVtxId);
1448     Clean(silseed_x);
1449     Clean(silseed_y);
1450     Clean(silseed_z);
1451     Clean(silseed_pt);
1452     Clean(silseed_eta);
1453     Clean(silseed_phi);
1454     Clean(silseed_eta_vtx);
1455     Clean(silseed_phi_vtx);
1456     Clean(silseed_crossing);
1457     Clean(silseed_charge);
1458     Clean(silseed_nMvtx);
1459     Clean(silseed_nIntt);
1460     Clean(silseed_clusterKeys);
1461     Clean(silseed_cluster_layer);
1462     Clean(silseed_cluster_globalX);
1463     Clean(silseed_cluster_globalY);
1464     Clean(silseed_cluster_globalZ);
1465     Clean(silseed_cluster_phi);
1466     Clean(silseed_cluster_eta);
1467     Clean(silseed_cluster_r);
1468     Clean(silseed_cluster_phiSize);
1469     Clean(silseed_cluster_zSize);
1470     Clean(silseed_cluster_strobeID);
1471     Clean(silseed_cluster_timeBucketID);
1472     Clean(silseed_ngmvtx);
1473     Clean(silseed_ngintt);
1474     Clean(silseed_cluster_gcluster_key);
1475     Clean(silseed_cluster_gcluster_layer);
1476     Clean(silseed_cluster_gcluster_X);
1477     Clean(silseed_cluster_gcluster_Y);
1478     Clean(silseed_cluster_gcluster_Z);
1479     Clean(silseed_cluster_gcluster_r);
1480     Clean(silseed_cluster_gcluster_phi);
1481     Clean(silseed_cluster_gcluster_eta);
1482     Clean(silseed_cluster_gcluster_edep);
1483     Clean(silseed_cluster_gcluster_adc);
1484     Clean(silseed_cluster_gcluster_phiSize);
1485     Clean(silseed_cluster_gcluster_zSize);
1486
1487     Clean(clusterKey);
1488     Clean(cluster_layer);
1489     Clean(cluster_chip);
1490     Clean(cluster_stave);
1491     Clean(cluster_globalX);
1492     Clean(cluster_globalY);
1493     Clean(cluster_globalZ);
1494     Clean(cluster_phi);
1495     Clean(cluster_eta);
1496     Clean(cluster_r);
1497     Clean(cluster_phiSize);
1498     Clean(cluster_zSize);
1499     Clean(cluster_adc);
1500     Clean(cluster_timeBucketID);
1501     Clean(cluster_ladderZId);
1502     Clean(cluster_ladderPhiId);
1503     Clean(cluster_LocalX);
1504     Clean(cluster_LocalY);
1505     Clean(cluster_matchedG4P_trackID);
1506     Clean(cluster_matchedG4P_PID);
1507     Clean(cluster_matchedG4P_E);
1508     Clean(cluster_matchedG4P_pT);
1509     Clean(cluster_matchedG4P_eta);
1510     Clean(cluster_matchedG4P_phi);
1511
1512     Clean(mvtx_seedcluster_key);
1513     Clean(mvtx_seedcluster_layer);
1514     Clean(mvtx_seedcluster_chip);
1515     Clean(mvtx_seedcluster_stave);
1516     Clean(mvtx_seedcluster_globalX);
1517     Clean(mvtx_seedcluster_globalY);
1518     Clean(mvtx_seedcluster_globalZ);
1519     Clean(mvtx_seedcluster_phi);
1520     Clean(mvtx_seedcluster_eta);
1521     Clean(mvtx_seedcluster_r);
1522     Clean(mvtx_seedcluster_phiSize);
1523     Clean(mvtx_seedcluster_zSize);
1524     Clean(mvtx_seedcluster_strobeID);
1525     Clean(mvtx_seedcluster_matchedcrossing);
1526     Clean(mvtx_seedcluster_hitX);
1527     Clean(mvtx_seedcluster_hitY);
1528     Clean(mvtx_seedcluster_hitZ);
1529     Clean(mvtx_seedcluster_hitrow);
1530     Clean(mvtx_seedcluster_hitcol);
1531
1532     Clean(mvtx_seedcluster_matchedG4P_trackID);
1533     Clean(mvtx_seedcluster_matchedG4P_PID);
1534     Clean(mvtx_seedcluster_matchedG4P_E);
1535     Clean(mvtx_seedcluster_matchedG4P_pT);
1536     Clean(mvtx_seedcluster_matchedG4P_eta);
1537     Clean(mvtx_seedcluster_matchedG4P_phi);
1538     Clean(mvtx_seedcluster_matchedG4P_ancestor_trackID);
1539     Clean(mvtx_seedcluster_matchedG4P_ancestor_PID);
1540
1541     Clean(PrimaryPHG4Ptcl_pT);
1542     Clean(PrimaryPHG4Ptcl_eta);
1543     Clean(PrimaryPHG4Ptcl_phi);
1544     Clean(PrimaryPHG4Ptcl_E);
1545     Clean(PrimaryPHG4Ptcl_PID);
1546     Clean(PrimaryPHG4Ptcl_trackID);
1547     Clean(PrimaryPHG4Ptcl_ParticleClass);
1548     Clean(PrimaryPHG4Ptcl_isStable);
1549     Clean(PrimaryPHG4Ptcl_charge);
1550     Clean(PrimaryPHG4Ptcl_isChargedHadron);
1551     Clean(PrimaryPHG4Ptcl_ancestor_trackID);
1552     Clean(PrimaryPHG4Ptcl_ancestor_PID);
1553     Clean(PrimaryPHG4Ptcl_truthcluster_X);
1554     Clean(PrimaryPHG4Ptcl_truthcluster_Y);
1555     Clean(PrimaryPHG4Ptcl_truthcluster_Z);
1556     Clean(PrimaryPHG4Ptcl_truthcluster_edep);
1557     Clean(PrimaryPHG4Ptcl_truthcluster_adc);
1558     Clean(PrimaryPHG4Ptcl_truthcluster_r);
1559     Clean(PrimaryPHG4Ptcl_truthcluster_phi);
1560     Clean(PrimaryPHG4Ptcl_truthcluster_eta);
1561     Clean(PrimaryPHG4Ptcl_truthcluster_phisize);
1562     Clean(PrimaryPHG4Ptcl_truthcluster_zsize);
1563     Clean(PrimaryPHG4Ptcl_recocluster_globalX);
1564     Clean(PrimaryPHG4Ptcl_recocluster_globalY);
1565     Clean(PrimaryPHG4Ptcl_recocluster_globalZ);
1566     Clean(PrimaryPHG4Ptcl_recocluster_r);
1567     Clean(PrimaryPHG4Ptcl_recocluster_phi);
1568     Clean(PrimaryPHG4Ptcl_recocluster_eta);
1569     Clean(PrimaryPHG4Ptcl_recocluster_phisize);
1570     Clean(PrimaryPHG4Ptcl_recocluster_zsize);
1571     Clean(PrimaryPHG4Ptcl_recocluster_adc);
1572
1573     Clean(sPHENIXPrimary_pT);
1574     Clean(sPHENIXPrimary_eta);
1575     Clean(sPHENIXPrimary_phi);
1576     Clean(sPHENIXPrimary_E);
1577     Clean(sPHENIXPrimary_PID);
1578     Clean(sPHENIXPrimary_trackID);
1579     Clean(sPHENIXPrimary_ParticleClass);
1580     Clean(sPHENIXPrimary_isStable);
1581     Clean(sPHENIXPrimary_charge);
1582     Clean(sPHENIXPrimary_isChargedHadron);
1583     Clean(sPHENIXPrimary_ancestor_trackID);
1584     Clean(sPHENIXPrimary_ancestor_PID);
1585     Clean(sPHENIXPrimary_truthcluster_X);
1586     Clean(sPHENIXPrimary_truthcluster_Y);
1587     Clean(sPHENIXPrimary_truthcluster_Z);
1588     Clean(sPHENIXPrimary_truthcluster_edep);
1589     Clean(sPHENIXPrimary_truthcluster_adc);
1590     Clean(sPHENIXPrimary_truthcluster_r);
1591     Clean(sPHENIXPrimary_truthcluster_phi);
1592     Clean(sPHENIXPrimary_truthcluster_eta);
1593     Clean(sPHENIXPrimary_truthcluster_phisize);
1594     Clean(sPHENIXPrimary_truthcluster_zsize);
1595     Clean(sPHENIXPrimary_recocluster_globalX);
1596     Clean(sPHENIXPrimary_recocluster_globalY);
1597     Clean(sPHENIXPrimary_recocluster_globalZ);
1598     Clean(sPHENIXPrimary_recocluster_r);
1599     Clean(sPHENIXPrimary_recocluster_phi);
1600     Clean(sPHENIXPrimary_recocluster_eta);
1601     Clean(sPHENIXPrimary_recocluster_phisize);
1602     Clean(sPHENIXPrimary_recocluster_zsize);
1603     Clean(sPHENIXPrimary_recocluster_adc);
1604
1605     Clean(AllPHG4Ptcl_pT);
1606     Clean(AllPHG4Ptcl_eta);
1607     Clean(AllPHG4Ptcl_phi);
1608     Clean(AllPHG4Ptcl_E);
1609     Clean(AllPHG4Ptcl_PID);
1610     Clean(AllPHG4Ptcl_trackID);
1611     Clean(AllPHG4Ptcl_ancestor_trackID);
1612     Clean(AllPHG4Ptcl_ancestor_PID);
1613     Clean(AllPHG4Ptcl_truthcluster_X);
1614     Clean(AllPHG4Ptcl_truthcluster_Y);
1615     Clean(AllPHG4Ptcl_truthcluster_Z);
1616     Clean(AllPHG4Ptcl_truthcluster_edep);
1617     Clean(AllPHG4Ptcl_truthcluster_adc);
1618     Clean(AllPHG4Ptcl_truthcluster_r);
1619     Clean(AllPHG4Ptcl_truthcluster_phi);
1620     Clean(AllPHG4Ptcl_truthcluster_eta);
1621     Clean(AllPHG4Ptcl_truthcluster_phisize);
1622     Clean(AllPHG4Ptcl_truthcluster_zsize);
1623     Clean(AllPHG4Ptcl_recocluster_globalX);
1624     Clean(AllPHG4Ptcl_recocluster_globalY);
1625     Clean(AllPHG4Ptcl_recocluster_globalZ);
1626     Clean(AllPHG4Ptcl_recocluster_r);
1627     Clean(AllPHG4Ptcl_recocluster_phi);
1628     Clean(AllPHG4Ptcl_recocluster_eta);
1629     Clean(AllPHG4Ptcl_recocluster_phisize);
1630     Clean(AllPHG4Ptcl_recocluster_zsize);
1631     Clean(AllPHG4Ptcl_recocluster_adc);
1632 }
1633
1634 //____________________________________________________________________________..
1635 int VertexCompare::ResetEvent(PHCompositeNode *topNode) { return Fun4AllReturnCodes::EVENT_OK; }
1636
1637 //____________________________________________________________________________..
1638 int VertexCompare::EndRun(const int runnumber) { return Fun4AllReturnCodes::EVENT_OK; }
1639
1640 //____________________________________________________________________________..
1641 int VertexCompare::End(PHCompositeNode *topNode)
1642 {
1643     outFile->Write("", TObject::kOverwrite);
1644     outFile->Close();
1645     delete outFile;
1646
1647     return Fun4AllReturnCodes::EVENT_OK;
1648 }
1649
1650 //____________________________________________________________________________..
1651 int VertexCompare::Reset(PHCompositeNode *topNode) { return Fun4AllReturnCodes::EVENT_OK; }
1652
1653 //____________________________________________________________________________..
1654 void VertexCompare::Print(const std::string &what) const {}