File indexing completed on 2026-08-31 08:21:25
0001 #ifndef TPCTRACKRECO_TPCCROSSINGFINDER_H
0002 #define TPCTRACKRECO_TPCCROSSINGFINDER_H
0003
0004
0005 #include <fun4all/SubsysReco.h>
0006
0007 #include "TpcCrossingDecision.h"
0008 #include <trackbase/TrkrDefs.h>
0009
0010 #include <array>
0011 #include <limits>
0012 #include <map>
0013 #include <set>
0014 #include <string>
0015 #include <utility>
0016 #include <vector>
0017
0018 class IdealPadMap;
0019 class PHCompositeNode;
0020 class PHGarfield;
0021 class PHG4TpcGeomContainer;
0022 class SvtxVertexMap;
0023 class Tpc_AssembledTrack;
0024 class Tpc_AssembledTrackContainer;
0025 class TpcCrossingDecisionContainerv1;
0026 class TrkrClusterContainer;
0027 class TrkrHitSetContainer;
0028
0029 class TpcCrossingFinder : public SubsysReco
0030 {
0031 public:
0032 explicit TpcCrossingFinder(const std::string& name = "TpcCrossingFinder");
0033 ~TpcCrossingFinder() override;
0034
0035 int InitRun(PHCompositeNode*) override;
0036 int process_event(PHCompositeNode*) override;
0037
0038 static constexpr unsigned int NPhiSamples = 3;
0039
0040 void setInputNodeName(const std::string& n) { m_inputNodeName = n; }
0041 void setOutputNodeName(const std::string& n) { m_outputNodeName = n; }
0042 void setVertexMapNodeName(const std::string& n) { m_vertexMapNodeName = n; }
0043 void setT0(double v) { m_t0 = v; }
0044 void setTpcAdcClock(double v) { m_tpcAdcClock = v; }
0045 void setCrossingPeriodNs(double v) { m_crossingPeriodNs = v; }
0046 void setReverseDriftStepNs(double v) { m_reverseDriftStepNs = v; }
0047 void setKEffSide0(double v) { m_kEffSide0 = v; }
0048 void setKEffSide1(double v) { m_kEffSide1 = v; }
0049 void setCMVoltageDefault(double v) { m_cmVoltageDefault = v; }
0050 void setRequireSiliconVertex(bool v) { m_requireSiliconVertex = v; }
0051 void setResolveAmbiguousWithoutVertex(bool v) { m_resolveAmbiguousWithoutVertex = v; }
0052 void setPreferTriggeredCrossing(bool v) { m_preferTriggeredCrossing = v; }
0053 void setTriggeredCrossing(short v) { m_triggeredCrossing = v; }
0054 void setCollisionZ(double v) { m_collisionZ = v; }
0055 void setMaxVertexDz(double v) { m_maxVertexDz = v; }
0056 void setMaxTier2BeamlineZ(double v) { m_maxTier2BeamlineZ = v; }
0057 void setMaxCandidateVertexZ(double v) { m_maxCandidateVertexZ = v; }
0058 void setMinBestSecondSeparation(double v) { m_minBestSecondSeparation = v; }
0059 void setTpcGeometryTolerance(double radial_cm, double z_cm, double central_membrane_cm)
0060 {
0061 m_radialTolerance = radial_cm;
0062 m_zTolerance = z_cm;
0063 m_centralMembraneTolerance = central_membrane_cm;
0064 }
0065 void setTpcHalfLength(double v) { m_tpcHalfLength = v; }
0066 void setUseSurveyGeometry(bool v) { use_survey_geometry = v; }
0067 void setMoveTpc(double x, double y, double z) { m_tpcMove = {{x, y, z}}; }
0068 void setRotateTpc(unsigned int index, double x, double y, double z)
0069 {
0070 if (index < m_tpcRotations.size()) m_tpcRotations[index] = {{x, y, z}};
0071 }
0072 void setStartZ(double south_z, double north_z)
0073 {
0074 m_startZSouth = south_z;
0075 m_startZNorth = north_z;
0076 }
0077
0078 private:
0079 struct DriftPoint
0080 {
0081 float delta_r {0.0F};
0082 float delta_phi {0.0F};
0083 float z {0.0F};
0084 };
0085
0086 struct DriftPolyline
0087 {
0088 double phi {0.0};
0089 std::vector<DriftPoint> points;
0090 };
0091
0092 struct Point
0093 {
0094 TrkrDefs::hitsetkey hitsetkey {0};
0095 TrkrDefs::hitkey hitkey {0};
0096 unsigned int layer {0};
0097 unsigned int side {0};
0098 unsigned int pad {0};
0099 unsigned int tbin {0};
0100 double x {0.0};
0101 double y {0.0};
0102 double z {0.0};
0103 };
0104
0105 struct SiliconVertexHypothesis
0106 {
0107 short crossing {0};
0108 unsigned int vertex_id {0};
0109 double x {0.0};
0110 double y {0.0};
0111 double z {0.0};
0112 double sigma_z {0.0};
0113 unsigned int ntracks {0};
0114 };
0115
0116 struct ZFitResult
0117 {
0118 bool valid {false};
0119 double slope {0.0};
0120 double intercept {0.0};
0121 double chi2 {0.0};
0122 int ndf {-1};
0123 double s_at_pca {0.0};
0124 double minimum_radius {0.0};
0125 double z_at_pca {0.0};
0126 double z_at_r0 {0.0};
0127 std::vector<Point> points;
0128 std::vector<float> path_length;
0129 };
0130
0131 struct Candidate
0132 {
0133 short crossing {0};
0134 bool tpc_valid {false};
0135 bool has_silicon_vertex {false};
0136 bool vertex_compatible {false};
0137 unsigned int silicon_vertex_id {0};
0138 double tpc_z0 {0.0};
0139 double silicon_vertex_z {0.0};
0140 double delta_z {0.0};
0141 unsigned char rejection_status {0};
0142 unsigned char confidence_tier {std::numeric_limits<unsigned char>::max()};
0143 double confidence_score {std::numeric_limits<double>::quiet_NaN()};
0144 TpcCrossingCandidate qa;
0145 };
0146
0147 int getNodes(PHCompositeNode*);
0148 int createNodes(PHCompositeNode*);
0149 void configure_garfield(PHGarfield* garfield) const;
0150 bool build_drift_lookup();
0151 bool sample_drift_lookup(unsigned int layer,
0152 unsigned int side,
0153 unsigned int pad,
0154 unsigned int tbin,
0155 short crossing,
0156 double& x,
0157 double& y,
0158 double& z) const;
0159 bool make_xyz_point(TrkrDefs::hitsetkey hsk, TrkrDefs::hitkey hk, short crossing, Point& p) const;
0160 bool find_time_extrema(const Tpc_AssembledTrack* track,
0161 TrkrDefs::hitsetkey& min_hsk,
0162 TrkrDefs::hitkey& min_hk,
0163 TrkrDefs::hitsetkey& max_hsk,
0164 TrkrDefs::hitkey& max_hk) const;
0165 std::set<short> get_available_crossings() const;
0166 std::set<short> get_intt_crossings() const;
0167 std::map<short, std::vector<SiliconVertexHypothesis>> get_vertices_by_crossing() const;
0168 std::vector<std::pair<TrkrDefs::hitsetkey, TrkrDefs::hitkey>> select_representatives(const Tpc_AssembledTrack* track,
0169 TrkrDefs::hitsetkey min_hsk,
0170 TrkrDefs::hitkey min_hk,
0171 TrkrDefs::hitsetkey max_hsk,
0172 TrkrDefs::hitkey max_hk) const;
0173 bool estimate_tpc_z0(std::vector<Point>& points, double& z0) const;
0174 ZFitResult estimate_tpc_z0_diagnostics(std::vector<Point> points) const;
0175 bool point_in_tpc(const Point& p) const;
0176 bool point_in_correct_side(const Point& p) const;
0177 Candidate test_candidate(const Tpc_AssembledTrack* track,
0178 short crossing,
0179 TrkrDefs::hitsetkey min_hsk,
0180 TrkrDefs::hitkey min_hk,
0181 TrkrDefs::hitsetkey max_hsk,
0182 TrkrDefs::hitkey max_hk,
0183 const std::map<short, std::vector<SiliconVertexHypothesis>>& vertices_by_crossing) const;
0184 static unsigned int drift_lookup_index(unsigned int layer_index, unsigned int side, unsigned int sector, unsigned int sample);
0185
0186 std::string m_inputNodeName {"TPC_ASSEMBLEDTRACKS"};
0187 std::string m_outputNodeName {"TPC_CROSSING_DECISIONS"};
0188 std::string m_vertexMapNodeName {"SvtxVertexMap"};
0189
0190 Tpc_AssembledTrackContainer* m_assembledTracks {nullptr};
0191 TpcCrossingDecisionContainerv1* m_decisions {nullptr};
0192 TrkrHitSetContainer* m_hits {nullptr};
0193 TrkrClusterContainer* m_clusterMap {nullptr};
0194 PHG4TpcGeomContainer* m_geomContainerTpc {nullptr};
0195 SvtxVertexMap* m_vertexMap {nullptr};
0196 IdealPadMap* m_idealPadMap {nullptr};
0197 PHGarfield* m_garfield {nullptr};
0198
0199 std::array<DriftPolyline, 48 * 2 * 12 * NPhiSamples> m_driftLookup;
0200 unsigned int m_event {0};
0201 double m_maxLookupTimeNs {0.0};
0202
0203 double m_t0 {8};
0204 double m_tpcAdcClock {56.881262};
0205 double m_crossingPeriodNs {106.56};
0206 double m_reverseDriftStepNs {56.881262};
0207 double m_startZSouth {-102.325};
0208 double m_startZNorth {102.325};
0209 double m_tpcHalfLength {105.5};
0210 double m_radialTolerance {2.0};
0211 double m_zTolerance {1.0};
0212 double m_centralMembraneTolerance {1.0};
0213 double m_cmVoltageDefault {380.0};
0214 double m_kEffSide0 {0.0};
0215 double m_kEffSide1 {-1.5};
0216 double m_maxVertexDz {2.0};
0217 double m_minBestSecondSeparation {0.4};
0218 double m_collisionZ {0.0};
0219 double m_maxTier2BeamlineZ {40.0};
0220 double m_maxCandidateVertexZ {20.0};
0221 bool m_requireSiliconVertex {false};
0222 bool m_resolveAmbiguousWithoutVertex {true};
0223 bool m_preferTriggeredCrossing {false};
0224 short m_triggeredCrossing {0};
0225 bool use_survey_geometry {false};
0226 std::array<double, 3> m_tpcMove {{0.0, 0.0, 0.0}};
0227 std::array<std::array<double, 3>, 2> m_tpcRotations {{{{0.0, 0.0, 0.0}}, {{0.0, 0.0, 0.0}}}};
0228 };
0229
0230 #endif