ParticleGen-Jingyu/physiTuto/TruthToSvtxTrack.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 TruthToSvtxTrack.h.
0007 //
0008 // TruthToSvtxTrack(const std::string &name = "TruthToSvtxTrack")
0009 // everything is keyed to TruthToSvtxTrack, 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 // TruthToSvtxTrack::~TruthToSvtxTrack()
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 TruthToSvtxTrack::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 TruthToSvtxTrack::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 TruthToSvtxTrack::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 TruthToSvtxTrack::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 TruthToSvtxTrack::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 TruthToSvtxTrack::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 TruthToSvtxTrack::Reset(PHCompositeNode *topNode)
0057 // not really used - it is called before the dtor is called
0058 //
0059 // void TruthToSvtxTrack::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 "TruthToSvtxTrack.h"
0065
0066 #include <globalvertex/GlobalVertexMap.h>
0067 #include <globalvertex/SvtxVertexMap.h>
0068 #include <ffarawobjects/Gl1Packet.h>
0069 #include <trackbase_historic/SvtxPHG4ParticleMap_v1.h>
0070 #include <trackbase_historic/SvtxTrack.h>
0071 #include <trackbase_historic/SvtxTrack_v1.h>
0072 #include <trackbase_historic/SvtxTrack_v2.h>
0073 #include <trackbase_historic/SvtxTrackMap.h>
0074 #include <trackbase_historic/SvtxTrackMap_v1.h>
0075 #include <trackbase_historic/SvtxTrackState_v1.h>
0076 #include <trackbase_historic/TrackSeed.h>
0077
0078 #include <calobase/RawClusterContainer.h>
0079 #include <calobase/RawTowerGeomContainer.h>
0080 #include <calobase/RawCluster.h>
0081 #include <calobase/RawClusterUtility.h>
0082 #include <calobase/RawTowerDefs.h>
0083 #include <calobase/RawTowerGeom.h>
0084 #include <calobase/RawTowerGeomv5.h>
0085 #include <calobase/TowerInfoContainer.h>
0086 #include <calobase/TowerInfo.h>
0087 #include <calobase/TowerInfoDefs.h>
0088
0089 #include <ffaobjects/EventHeaderv1.h>
0090
0091 #include <fun4all/Fun4AllReturnCodes.h>
0092
0093 #include <globalvertex/GlobalVertex.h>
0094 #include <globalvertex/GlobalVertexMap.h>
0095 #include <globalvertex/MbdVertexMap.h>
0096 #include <globalvertex/MbdVertex.h>
0097 #include <globalvertex/SvtxVertexMap.h>
0098 #include <globalvertex/SvtxVertex.h>
0099
0100 #include <phool/getClass.h>
0101 #include <phool/PHCompositeNode.h>
0102
0103 #include <Acts/Geometry/GeometryIdentifier.hpp>
0104 #include <Acts/MagneticField/ConstantBField.hpp>
0105 #include <Acts/MagneticField/MagneticFieldProvider.hpp>
0106 #include <Acts/Surfaces/CylinderSurface.hpp>
0107 #include <Acts/Surfaces/PerigeeSurface.hpp>
0108 #include <Acts/Geometry/TrackingGeometry.hpp>
0109
0110 #include <CLHEP/Vector/ThreeVector.h>
0111 #include <math.h>
0112 #include <vector>
0113
0114 #include <TFile.h>
0115 #include <TTree.h>
0116 #include <TH1F.h>
0117 #include <TH2F.h>
0118 #include <TLorentzVector.h>
0119
0120 //____________________________________________________________________________..
0121 TruthToSvtxTrack::TruthToSvtxTrack(
0122     const std::string & name_in):
0123     SubsysReco(name_in)
0124 {
0125   std::cout << "TruthToSvtxTrack::TruthToSvtxTrack(const std::string &name) Calling ctor" << std::endl;
0126 }
0127
0128 //____________________________________________________________________________..
0129 TruthToSvtxTrack::~TruthToSvtxTrack()
0130 {
0131   std::cout << "TruthToSvtxTrack::~TruthToSvtxTrack() Calling dtor" << std::endl;
0132 }
0133
0134 //____________________________________________________________________________..
0135 int TruthToSvtxTrack::Init(PHCompositeNode *topNode)
0136 {
0137     std::cout << topNode << std::endl;
0138     std::cout << "TruthToSvtxTrack::Init(PHCompositeNode *topNode) Initializing" << std::endl;
0139
0140     return Fun4AllReturnCodes::EVENT_OK;
0141 }
0142
0143 //____________________________________________________________________________..
0144 int TruthToSvtxTrack::InitRun(PHCompositeNode * /*topNode*/)
0145 {
0146   std::cout << "TruthToSvtxTrack::InitRun(PHCompositeNode *topNode) Initializing for Run XXX" << std::endl;
0147   return Fun4AllReturnCodes::EVENT_OK;
0148 }
0149
0150 //____________________________________________________________________________..
0151 int TruthToSvtxTrack::process_event(PHCompositeNode * topNode) {
0152     createTracksFromTruth(topNode);
0153
0154     return Fun4AllReturnCodes::EVENT_OK;
0155 }
0156
0157 //____________________________________________________________________________..
0158 int TruthToSvtxTrack::ResetEvent(PHCompositeNode *topNode)
0159 {
0160     std::cout << "TruthToSvtxTrack::ResetEvent(PHCompositeNode *topNode) Resetting internal structures, prepare for next event" << std::endl;
0161
0162     return Fun4AllReturnCodes::EVENT_OK;
0163 }
0164
0165 //____________________________________________________________________________..
0166 int TruthToSvtxTrack::EndRun(const int runnumber)
0167 {
0168     std::cout << "TruthToSvtxTrack::EndRun(const int runnumber) Ending Run for Run " << runnumber << std::endl;
0169
0170     return Fun4AllReturnCodes::EVENT_OK;
0171 }
0172
0173 //____________________________________________________________________________..
0174 int TruthToSvtxTrack::End(PHCompositeNode *topNode)
0175 {
0176   std::cout << "TruthToSvtxTrack::End(PHCompositeNode *topNode) This is the End..." << std::endl;
0177
0178   return Fun4AllReturnCodes::EVENT_OK;
0179 }
0180
0181 //____________________________________________________________________________..
0182 int TruthToSvtxTrack::Reset(PHCompositeNode *topNode)
0183 {
0184  std::cout << "TruthToSvtxTrack::Reset(PHCompositeNode *topNode) being Reset" << std::endl;
0185   return Fun4AllReturnCodes::EVENT_OK;
0186 }
0187
0188 //____________________________________________________________________________..
0189 void TruthToSvtxTrack::Print(const std::string &what) const
0190 {
0191   std::cout << "TruthToSvtxTrack::Print(const std::string &what) const Printing info for " << what << std::endl;
0192 }
0193
0194 int TruthToSvtxTrack::createTracksFromTruth(PHCompositeNode* topNode)
0195 {
0196     m_truth_info = findNode::getClass<PHG4TruthInfoContainer>(topNode, "G4TruthInfo");
0197
0198     if (!m_truth_info)
0199     {
0200         std::cout << PHWHERE << "Error, can't find G4TruthInfo" << std::endl;
0201         exit(1);
0202     }
0203
0204     // 获取或创建 SvtxTrackMap 节点
0205     SvtxTrackMap* trackmap = findNode::getClass<SvtxTrackMap>(topNode, "SvtxTrackMap"); // "TruthTrackMap"
0206     if (!trackmap)
0207     {
0208         trackmap = new SvtxTrackMap_v1();
0209         PHNodeIterator iter(topNode);
0210         PHCompositeNode* dstNode = dynamic_cast<PHCompositeNode*>(iter.findFirst("PHCompositeNode", "DST"));
0211
0212         PHIODataNode<PHObject> *truthtracknode = new PHIODataNode<PHObject>(trackmap, "SvtxTrackMap", "PHObject");
0213         dstNode->addNode(truthtracknode);
0214     }
0215
0216     const auto prange = m_truth_info->GetPrimaryParticleRange();
0217     for (auto iter = prange.first; iter != prange.second; ++iter)
0218     {
0219         PHG4Particle* ptcl = iter->second;
0220         if (!ptcl) continue;
0221
0222         TLorentzVector p;
0223         p.SetPxPyPzE(ptcl->get_px(), ptcl->get_py(), ptcl->get_pz(), ptcl->get_e());
0224
0225         SvtxTrack_v2* track = new SvtxTrack_v2();
0226         track->set_id(trackmap->size());
0227         track->set_px(ptcl->get_px());
0228         track->set_py(ptcl->get_py());
0229         track->set_pz(ptcl->get_pz());
0230         track->set_charge(static_cast<int>(TDatabasePDG::Instance()->GetParticle(ptcl->get_pid())->Charge()));
0231         track->set_chisq(1.0);
0232         track->set_ndf(1);
0233         track->set_vertex_id(0);  // 可改为对应 vtx id
0234
0235         trackmap->insert(track);
0236     }
0237
0238     return Fun4AllReturnCodes::EVENT_OK;
0239 }
0240