packages/tpctrackreco/Tpc_AssembledTrackReco.h

0001 // Tell emacs that this is a C++ source
0002 //  -*- C++ -*-.
0003 #ifndef TPCTRACKRECO_TPCASSEMBLEDTRACKRECO_H
0004 #define TPCTRACKRECO_TPCASSEMBLEDTRACKRECO_H
0005
0006 #include <fun4all/SubsysReco.h>
0007 #include <trackbase/TrkrDefs.h>
0008
0009 class IdealPadMap;
0010
0011 #include <mutex>
0012 #include <string>
0013 #include <vector>
0014
0015 class PHCompositeNode;
0016 class TFile;
0017 class TTree;
0018 class TH1D;
0019 class TH2D;
0020
0021 class Tpc_ModuleTrackContainer;
0022 class Tpc_AssembledTrackContainer;
0023 class TrkrHitSetContainer;
0024
0025 class Tpc_AssembledTrackReco : public SubsysReco
0026 {
0027  public:
0028   explicit Tpc_AssembledTrackReco(const std::string& name = "Tpc_AssembledTrackReco",
0029                                   const std::string& filename = "Tpc_AssembledTracks.root");
0030   ~Tpc_AssembledTrackReco() override;
0031
0032   int Init(PHCompositeNode*) override;
0033   int InitRun(PHCompositeNode*) override;
0034   int process_event(PHCompositeNode*) override;
0035   int End(PHCompositeNode*) override;
0036
0037   void setInputNodeName(const std::string& n) { m_inputNodeName = n; }
0038   void setOutputNodeName(const std::string& n) { m_outputNodeName = n; }
0039   void setDebugOutputFileName(const std::string& n) { m_debugOutputFileName = n; }
0040
0041   void setConnectMaxLayerGap(unsigned int n) { m_connectMaxLayerGap = n; }
0042
0043   // Connection window in global phi radians and tbin units at the match radius.
0044   void setConnectWindow(double dphi, double dtbin)
0045   {
0046     m_connect_dphi = dphi;
0047     m_connect_dtbin = dtbin;
0048   }
0049
0050   // Slope windows: d(phi)/d(radius) and d(tbin)/d(radius).
0051   void setConnectSlopeWindow(double dphi_slope, double dtbin_slope)
0052   {
0053     m_connect_dphi_slope = dphi_slope;
0054     m_connect_dtbin_slope = dtbin_slope;
0055   }
0056
0057   void setUseSagittaPhiFit(bool v) { m_useSagittaPhiFit = v; }
0058
0059   void setSeedCovarianceDiagonal(double sigmaX, double sigmaY, double sigmaZ,
0060                                  double sigmaPx, double sigmaPy, double sigmaPz)
0061   {
0062     m_seedSigmaX = sigmaX;
0063     m_seedSigmaY = sigmaY;
0064     m_seedSigmaZ = sigmaZ;
0065     m_seedSigmaPx = sigmaPx;
0066     m_seedSigmaPy = sigmaPy;
0067     m_seedSigmaPz = sigmaPz;
0068   }
0069
0070  public:
0071   // ---------------------------------------------------------------
0072   // Piece: one Tpc_ModuleTrack refit from stored hit indices for temporary matching.
0073   // Fits are in ideal (radius, global phi) and (radius, timebin) coordinates.
0074   // ---------------------------------------------------------------
0075   struct Piece
0076   {
0077     Piece();
0078
0079     unsigned int source_index;
0080     unsigned int source_track_id;
0081     unsigned int event;
0082     unsigned int region;
0083     unsigned int sector;
0084     int side;
0085
0086     unsigned int first_layer;
0087     unsigned int last_layer;
0088     unsigned int nblobs;
0089     unsigned int nrawhits;
0090
0091     double phi_slope;
0092     double phi_intercept;
0093     double phi_S;
0094     double phi_x0;
0095     double phi_invR;
0096     double phi_theta;
0097     double phi_bline;
0098     bool phi_sagitta_ok;
0099
0100     double tbin_slope;
0101     double tbin_intercept;
0102
0103     std::vector<double> radius_values;
0104     std::vector<double> phi_values;
0105     std::vector<double> tbin_values;
0106     std::vector<double> weights;
0107
0108     std::vector<TrkrDefs::hitsetkey> hitsetkeys;
0109     std::vector<TrkrDefs::hitkey> hitkeys;
0110   };
0111
0112   // ---------------------------------------------------------------
0113   // Candidate: a growing assembled track built from one or more Pieces.
0114   // Fits are temporary and are not copied into the output Tpc_AssembledTrack object.
0115   // ---------------------------------------------------------------
0116   struct SeedParameters
0117   {
0118     bool ok{false};
0119     double x{0.0};
0120     double y{0.0};
0121     double z{0.0};
0122     double px{0.0};
0123     double py{0.0};
0124     double pz{0.0};
0125     double cov[6][6]{};
0126   };
0127
0128   struct Candidate
0129   {
0130     Candidate();
0131
0132     unsigned int event;
0133     int side;
0134     unsigned int first_layer;
0135     unsigned int last_layer;
0136     unsigned int first_sector;
0137     unsigned int last_sector;
0138     unsigned int first_region;
0139     unsigned int last_region;
0140     unsigned int nsegments;
0141     unsigned int nblobs;
0142     unsigned int nrawhits;
0143
0144     double phi_slope;
0145     double phi_intercept;
0146     double phi_S;
0147     double phi_x0;
0148     double phi_invR;
0149     double phi_theta;
0150     double phi_bline;
0151     bool phi_sagitta_ok;
0152
0153     double tbin_slope_r;
0154     double tbin_intercept_r;
0155     double chi2_phi;
0156     double chi2_tbin;
0157     int ndof_phi;
0158     int ndof_tbin;
0159
0160     std::vector<unsigned int> piece_indices;
0161     std::vector<TrkrDefs::hitsetkey> hitsetkeys;
0162     std::vector<TrkrDefs::hitkey> hitkeys;
0163   };
0164
0165  private:
0166   int getNodes(PHCompositeNode*);
0167   int createNodes(PHCompositeNode*);
0168   void reset_tree_vars();
0169   void create_debug_histograms();
0170   void write_debug_histograms();
0171
0172   bool make_piece(unsigned int source_index, Piece& p) const;
0173
0174   double predict_phi(const Piece& p, double radius) const;
0175   double predict_phi(const Candidate& c, double radius) const;
0176   double predict_phi_slope(const Piece& p, double radius) const;
0177   double predict_phi_slope(const Candidate& c, double radius) const;
0178
0179   // Refit a growing candidate for temporary matching only.
0180   bool refit_candidate(const std::vector<Piece>& pieces,
0181                        const std::vector<unsigned int>& piece_indices,
0182                        Candidate& c) const;
0183
0184   // Check whether piece b can be appended to candidate a.
0185   bool candidates_can_connect(const Candidate& a,
0186                               const Piece& b,
0187                               double& score,
0188                               double& b_phi_intercept_shifted) const;
0189   bool candidates_can_connect(const Candidate& a,
0190                               const Candidate& b,
0191                               double& score) const;
0192
0193   void connect_sector_pieces(const std::vector<Piece>& pieces,
0194                              int side,
0195                              unsigned int sector,
0196                              std::vector<Candidate>& output) const;
0197   void connect_side_candidates(const std::vector<Piece>& pieces,
0198                                const std::vector<Candidate>& seeds,
0199                                int side,
0200                                std::vector<Candidate>& output) const;
0201   SeedParameters make_seed_parameters(const Candidate& c) const;
0202
0203   // ---------------------------------------------------------------
0204   std::string m_outputFileName;
0205   std::string m_debugOutputFileName;
0206   std::string m_inputNodeName;
0207   std::string m_outputNodeName;
0208
0209   TFile* m_outputFile;
0210   TFile* m_debugOutputFile;
0211   TTree* m_tree;
0212
0213   Tpc_ModuleTrackContainer* m_tpcModuleTrackContainer;
0214   Tpc_AssembledTrackContainer* m_assembledTrackContainer;
0215   TrkrHitSetContainer* m_hits;
0216
0217   int m_event;
0218   IdealPadMap* m_idealPadMap;
0219
0220   unsigned int m_connectMaxLayerGap;
0221   double m_connect_dphi;
0222   double m_connect_dtbin;
0223   double m_connect_dphi_slope;
0224   double m_connect_dtbin_slope;
0225   bool m_useSagittaPhiFit;
0226   double m_seedSigmaX;
0227   double m_seedSigmaY;
0228   double m_seedSigmaZ;
0229   double m_seedSigmaPx;
0230   double m_seedSigmaPy;
0231   double m_seedSigmaPz;
0232
0233   // Debug matching histograms
0234   TH1D* m_h_dphi;
0235   TH1D* m_h_dtbin;
0236   TH1D* m_h_dmphi;
0237   TH1D* m_h_dmtbin;
0238   TH1D* m_h_score;
0239   TH2D* m_h_dphi_vs_dtbin;
0240   TH2D* m_h_dmphi_vs_dmtbin;
0241   TH2D* m_h_dphi_vs_dmphi;
0242   TH2D* m_h_tbin_slope_vs_first_tbin;
0243   TH2D* m_h_tbin_slope_vs_last_tbin;
0244   TH2D* m_h_track_tbin_slope_vs_tbin_span_3modules;
0245   TH2D* m_h_track_tbin_slope_vs_first_tbin_3modules;
0246   TH2D* m_h_track_tbin_slope_vs_last_tbin_3modules;
0247   TH1D* m_h_layer_gap;
0248   TH1D* m_h_nsegments;
0249   TH1D* m_h_matched_sector_delta;
0250
0251   mutable std::mutex m_debugMutex;
0252
0253   // ---------------------------------------------------------------
0254   // TTree branches
0255   // ---------------------------------------------------------------
0256   int m_tree_event{};
0257   std::vector<unsigned int> m_tree_track_id;
0258   std::vector<int> m_tree_side;
0259   std::vector<unsigned int> m_tree_nsegments;
0260   std::vector<unsigned int> m_tree_nblobs;
0261   std::vector<unsigned int> m_tree_nrawhits;
0262   std::vector<unsigned int> m_tree_first_layer;
0263   std::vector<unsigned int> m_tree_last_layer;
0264   std::vector<unsigned int> m_tree_first_sector;
0265   std::vector<unsigned int> m_tree_last_sector;
0266   std::vector<unsigned int> m_tree_first_region;
0267   std::vector<unsigned int> m_tree_last_region;
0268
0269   std::vector<unsigned int> m_tree_source_assembled_track_id;
0270   std::vector<unsigned int> m_tree_source_inmodule_track_id;
0271   std::vector<unsigned int> m_tree_source_region;
0272   std::vector<unsigned int> m_tree_source_sector;
0273   std::vector<int> m_tree_source_side;
0274
0275   std::vector<unsigned int> m_tree_hit_assembled_track_id;
0276   std::vector<unsigned long long> m_tree_hit_hitsetkey;
0277   std::vector<unsigned long long> m_tree_hit_hitkey;
0278 };
0279 #endif