Warning, /coresoftware/offline/packages/trackreco/PHTpcTrackSeedCircleFit.cc.outdated is written in an unsupported language. File is not indexed.
0001 #include "PHTpcTrackSeedCircleFit.h"
0002
0003 #include <trackbase/TrkrDefs.h> // for cluskey, getTrkrId, tpcId
0004 #include <trackbase/TrkrClusterContainer.h>
0005 #include <trackbase/TrkrCluster.h>
0006 #include <trackbase_historic/TrackSeed.h>
0007 #include <trackbase_historic/TrackSeedContainer.h>
0008 #include <trackbase_historic/ActsTransformations.h>
0009
0010 #include <tpc/TpcDistortionCorrectionContainer.h>
0011
0012 #include <g4main/PHG4Hit.h> // for PHG4Hit
0013 #include <g4main/PHG4Particle.h> // for PHG4Particle
0014 #include <g4main/PHG4HitDefs.h> // for keytype
0015
0016 #include <fun4all/Fun4AllReturnCodes.h>
0017
0018 #include <phool/getClass.h>
0019 #include <phool/phool.h>
0020
0021 #if __cplusplus < 201402L
0022 #include <boost/make_unique.hpp>
0023 #endif
0024
0025 #include <TF1.h>
0026
0027 #include <climits> // for UINT_MAX
0028 #include <iostream> // for operator<<, basic_ostream
0029 #include <cmath> // for fabs, std::sqrt
0030 #include <set> // for _Rb_tree_const_iterator
0031 #include <utility> // for pair
0032 #include <memory>
0033
0034 namespace
0035 {
0036
0037 // square
0038 template<class T> inline constexpr T square( const T& x ) { return x*x; }
0039
0040 // radius
0041 template<class T> T get_r( const T& x, const T& y ) { return std::sqrt( square(x) + square(y) ); }
0042
0043 }
0044
0045 //____________________________________________________________________________..
0046 PHTpcTrackSeedCircleFit::PHTpcTrackSeedCircleFit(const std::string &name):
0047 SubsysReco(name)
0048 {}
0049
0050 //____________________________________________________________________________..
0051 int PHTpcTrackSeedCircleFit::InitRun(PHCompositeNode *topnode)
0052 {
0053 // get relevant nodes
0054 int ret = GetNodes( topnode );
0055 if( ret != Fun4AllReturnCodes::EVENT_OK ) return ret;
0056
0057 return Fun4AllReturnCodes::EVENT_OK; }
0058
0059 //____________________________________________________________________________..
0060 int PHTpcTrackSeedCircleFit::process_event(PHCompositeNode*)
0061 {
0062
0063 // _track_map contains the TPC seed track stubs
0064 // We want to associate these TPC track seeds with a collision vertex
0065 // All we need is to project the TPC clusters in Z to the beam line.
0066 // The TPC track seeds are given a position that is the PCA of the line and circle fit to the beam line
0067
0068 if(Verbosity() > 0)
0069 std::cout << PHWHERE << " TPC track map size " << _track_map->size() << std::endl;
0070
0071 unsigned int track_key = 0;
0072 // loop over the TPC track seeds
0073 for (auto phtrk_iter = _track_map->begin();
0074 phtrk_iter != _track_map->end();
0075 ++phtrk_iter)
0076 {
0077 TrackSeed* tracklet_tpc = *phtrk_iter;
0078
0079 if (Verbosity() > 1)
0080 {
0081 std::cout
0082 << __LINE__
0083 << ": Processing seed itrack: " << track_key
0084 << ": nhits: " << tracklet_tpc-> size_cluster_keys()
0085 << ": pT: " << tracklet_tpc->get_pt()
0086 << ": phi: " << tracklet_tpc->get_phi(_cluster_map, _surfmaps, _tGeometry)
0087 << ": eta: " << tracklet_tpc->get_eta()
0088 << std::endl;
0089 }
0090
0091 if(tracklet_tpc->size_cluster_keys() < 3)
0092 {
0093 if(Verbosity() > 3) std::cout << PHWHERE << " -- skip this tpc tracklet, not enough TPC clusters " << std::endl;
0094 continue; // skip to the next TPC tracklet
0095 }
0096
0097 /// Start at layer 7
0098 tracklet_tpc->lineFit(_cluster_map, _surfmaps, _tGeometry, 7, 58);
0099 tracklet_tpc->circleFitByTaubin(_cluster_map, _surfmaps, _tGeometry, 7, 58);
0100
0101 if(Verbosity() > 5)
0102 {
0103 std::cout << " new mom " << tracklet_tpc->get_p() << " new eta " << tracklet_tpc->get_eta()
0104 << " new phi " << tracklet_tpc->get_phi(_cluster_map, _surfmaps, _tGeometry) * 180.0 / M_PI << std::endl;
0105 }
0106
0107 track_key++;
0108 } // end loop over TPC track seeds
0109
0110 if(Verbosity() > 0)
0111 std::cout << " Final track map size " << _track_map->size() << std::endl;
0112
0113 if (Verbosity() > 0)
0114 std::cout << "PHTpcTrackSeedCircleFit::process_event(PHCompositeNode *topNode) Leaving process_event" << std::endl;
0115
0116 return Fun4AllReturnCodes::EVENT_OK;
0117 }
0118
0119 int PHTpcTrackSeedCircleFit::End(PHCompositeNode*)
0120 { return Fun4AllReturnCodes::EVENT_OK; }
0121
0122 int PHTpcTrackSeedCircleFit::GetNodes(PHCompositeNode* topNode)
0123 {
0124 _surfmaps = findNode::getClass<ActsSurfaceMaps>(topNode,"ActsSurfaceMaps");
0125 if(!_surfmaps)
0126 {
0127 std::cout << PHWHERE << "Error, can't find acts surface maps" << std::endl;
0128 return Fun4AllReturnCodes::ABORTEVENT;
0129 }
0130
0131 _tGeometry = findNode::getClass<ActsTrackingGeometry>(topNode,"ActsTrackingGeometry");
0132 if(!_tGeometry)
0133 {
0134 std::cout << PHWHERE << "Error, can't find acts tracking geometry" << std::endl;
0135 return Fun4AllReturnCodes::ABORTEVENT;
0136 }
0137
0138 if(_use_truth_clusters)
0139 _cluster_map = findNode::getClass<TrkrClusterContainer>(topNode, "TRKR_CLUSTER_TRUTH");
0140 else
0141 {
0142 _cluster_map = findNode::getClass<TrkrClusterContainer>(topNode, "TRKR_CLUSTER");
0143 }
0144 if (!_cluster_map)
0145 {
0146 std::cerr << PHWHERE << " ERROR: Can't find node TRKR_CLUSTER" << std::endl;
0147 return Fun4AllReturnCodes::ABORTEVENT;
0148 }
0149
0150 // tpc distortion correction
0151 _dcc = findNode::getClass<TpcDistortionCorrectionContainer>(topNode,"TpcDistortionCorrectionContainerStatic");
0152 if( _dcc )
0153 { std::cout << "PHTpcTrackSeedCircleFit::get_Nodes - found static TPC distortion correction container" << std::endl; }
0154
0155 _track_map = findNode::getClass<TrackSeedContainer>(topNode, _track_map_name);
0156 if (!_track_map)
0157 {
0158 std::cerr << PHWHERE << " ERROR: Can't find TrackSeedContainer" << std::endl;
0159 return Fun4AllReturnCodes::ABORTEVENT;
0160 }
0161
0162 return Fun4AllReturnCodes::EVENT_OK;
0163 }
0164
0165
0166 Acts::Vector3 PHTpcTrackSeedCircleFit::getGlobalPosition( TrkrDefs::cluskey key, TrkrCluster* cluster ) const
0167 {
0168 // get global position from Acts transform
0169 auto globalpos = _surfmaps->getGlobalPosition(key, cluster, _tGeometry);
0170
0171 // ADF: in streaming mode we need to add a step here to take care of the fact that we do not know the crossing yet
0172 // possibly move the track to point at z=0 to make distortion corrections (circularize the track) then move it back after the fit?
0173
0174 // check if TPC distortion correction are in place and apply
0175 if(_dcc) { globalpos = _distortionCorrection.get_corrected_position( globalpos, _dcc ); }
0176
0177 return globalpos;
0178 }
0179
