Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:17:16

0001 #ifndef TRACKRECO_PHACTSSILICONSEEDING_H
0002 #define TRACKRECO_PHACTSSILICONSEEDING_H
0003 
0004 #include <fun4all/SubsysReco.h>
0005 #include <trackbase/ActsGeometry.h>
0006 #include <trackbase/ClusterErrorPara.h>
0007 #include <trackbase/TrkrDefs.h>
0008 
0009 #include <Acts/Definitions/Algebra.hpp>
0010 #include <Acts/Definitions/Units.hpp>
0011 #include <Acts/Geometry/GeometryIdentifier.hpp>
0012 #include <Acts/Seeding/SeedFinder.hpp>
0013 
0014 #include <Acts/Utilities/GridBinFinder.hpp>
0015 #include <Acts/Seeding/SpacePointGrid.hpp>
0016 
0017 #include <trackbase/SpacePoint.h>
0018 
0019 #include <TFile.h>
0020 #include <TH1.h>
0021 #include <TH2.h>
0022 #include <TTree.h>
0023 #include <map>
0024 #include <string>
0025 
0026 class PHCompositeNode;
0027 class PHG4CylinderGeomContainer;
0028 class TrackSeed;
0029 class TrackSeedContainer;
0030 class TrkrCluster;
0031 class TrkrClusterContainer;
0032 class TrkrClusterIterationMapv1;
0033 class TrkrClusterCrossingAssoc;
0034 
0035 using GridSeeds = std::vector<std::vector<Acts::Seed<SpacePoint>>>;
0036 
0037 /**
0038  * This class runs the Acts seeder over the MVTX measurements
0039  * to create track stubs for the rest of the stub matching pattern
0040  * recognition. The module also projects the MVTX stubs to the INTT
0041  * to find possible matches in the INTT to the MVTX triplet.
0042  */
0043 class PHActsSiliconSeeding : public SubsysReco
0044 {
0045  public:
0046   PHActsSiliconSeeding(const std::string &name = "PHActsSiliconSeeding");
0047   ~PHActsSiliconSeeding() override;
0048   int Init(PHCompositeNode *topNode) override;
0049   int InitRun(PHCompositeNode *topNode) override;
0050   int process_event(PHCompositeNode *topNode) override;
0051   int End(PHCompositeNode *topNode) override;
0052 
0053   void setStrobeRange(const int low, const int high)
0054   {
0055     m_lowStrobeIndex = low;
0056     m_highStrobeIndex = high;
0057   }
0058   void setunc(float unc) { m_uncfactor = unc; }
0059   /// Set seeding with truth clusters
0060   void useTruthClusters(bool useTruthClusters)
0061   {
0062     m_useTruthClusters = useTruthClusters;
0063   }
0064 
0065   /// Output some diagnostic histograms
0066   void seedAnalysis(bool seedAnalysis)
0067   {
0068     m_seedAnalysis = seedAnalysis;
0069   }
0070 
0071   void setinttRPhiSearchWindow(const float win)
0072   {
0073     m_inttrPhiSearchWin = win;
0074   }
0075   void setinttZSearchWindow(const float &win)
0076   {
0077     m_inttzSearchWin = win;
0078   }
0079   void setmvtxRPhiSearchWindow(const float win)
0080   {
0081     m_mvtxrPhiSearchWin = win;
0082   }
0083   void setmvtxZSearchWindow(const float &win)
0084   {
0085     m_mvtxzSearchWin = win;
0086   }
0087   /// For each MVTX+INTT seed, take the best INTT hits and form
0088   /// 1 silicon seed per MVTX seed
0089   void cleanSeeds(bool cleanSeeds)
0090   {
0091     m_cleanSeeds = cleanSeeds;
0092   }
0093 
0094   void rMax(const float rMax)
0095   {
0096     m_rMax = rMax;
0097   }
0098   void rMin(const float rMin)
0099   {
0100     m_rMin = rMin;
0101   }
0102   void zMax(const float zMax)
0103   {
0104     m_zMax = zMax;
0105   }
0106   void zMin(const float zMin)
0107   {
0108     m_zMin = zMin;
0109   }
0110   void deltaRMax(const float deltaRMax)
0111   {
0112     m_deltaRMax = deltaRMax;
0113   }
0114   void cotThetaMax(const float cotThetaMax)
0115   {
0116     m_cotThetaMax = cotThetaMax;
0117   }
0118   void gridFactor(const float gridFactor)
0119   {
0120     m_gridFactor = gridFactor;
0121   }
0122   void sigmaScattering(const float sigma)
0123   {
0124     m_sigmaScattering = sigma;
0125   }
0126   void maxPtScattering(const float pt)
0127   {
0128     m_maxPtScattering = pt;
0129   }
0130   void sigmaError(const float sigma)
0131   {
0132     m_sigmaError = sigma;
0133   }
0134   void zalign(const float z)
0135   {
0136     m_zalign = z;
0137   }
0138   void ralign(const float r)
0139   {
0140     m_ralign = r;
0141   }
0142   void tolerance(const float tolerance)
0143   {
0144     m_tolerance = tolerance;
0145   }
0146   void helixcut(const float cut)
0147   {
0148     m_helixcut = cut;
0149   }
0150   void bfield(const float field)
0151   {
0152     m_bField = field;
0153   }
0154   void minpt(const float pt)
0155   {
0156     m_minSeedPt = pt;
0157   }
0158 
0159   /// A function to run the seeder with large (true)
0160   /// or small (false) grid spacing
0161   void largeGridSpacing(const bool spacing);
0162 
0163   void set_track_map_name(const std::string &map_name) { _track_map_name = map_name; }
0164   void iteration(int iter) { m_nIteration = iter; }
0165   void searchInIntt() { m_searchInIntt = true; }
0166 
0167  private:
0168   int getNodes(PHCompositeNode *topNode);
0169   int createNodes(PHCompositeNode *topNode);
0170 
0171   void runSeeder();
0172 
0173   /// Configure the seeding parameters for Acts. There
0174   /// are a number of tunable parameters for the seeder here
0175   void configureSeeder();
0176   void configureSPGrid();
0177   Acts::SeedFilterConfig configureSeedFilter();
0178 
0179   /// Take final seeds and fill the TrackSeedContainer
0180   void makeSvtxTracks(const GridSeeds &seedVector);
0181 
0182   /// Create a seeding space point out of an Acts::SourceLink
0183   SpacePointPtr makeSpacePoint(
0184       const Surface &surf,
0185       const TrkrDefs::cluskey,
0186       //    const TrkrCluster* clus);
0187       TrkrCluster *clus);
0188 
0189   /// Get all space points for the seeder
0190   std::vector<const SpacePoint *> getSiliconSpacePoints(Acts::Extent &rRangeSPExtent,
0191                                                         const int strobe);
0192   void printSeedConfigs(Acts::SeedFilterConfig &sfconfig);
0193 
0194   /// Projects circle fit to radii to find possible MVTX/INTT clusters
0195   /// belonging to track stub
0196   std::vector<TrkrDefs::cluskey> findMatches(
0197       std::vector<Acts::Vector3> &clusters,
0198       std::vector<TrkrDefs::cluskey> &keys,
0199       TrackSeed &seed);
0200 
0201   std::vector<TrkrDefs::cluskey> matchInttClusters(std::vector<Acts::Vector3> &clusters,
0202                                                    TrackSeed &seed,
0203                                                    const double xProj[],
0204                                                    const double yProj[],
0205                                                    const double zProj[]);
0206   short int getCrossingIntt(TrackSeed &si_track);
0207   std::vector<short int> getInttCrossings(TrackSeed &si_track);
0208 
0209   void createHistograms();
0210   void writeHistograms();
0211   double normPhi2Pi(const double phi);
0212   void clearTreeVariables();
0213 
0214   TrkrClusterCrossingAssoc *_cluster_crossing_map = nullptr;
0215   TTree *m_tree = nullptr;
0216   int m_seedid = std::numeric_limits<int>::quiet_NaN();
0217   std::vector<float> m_mvtxgx = {};
0218   std::vector<float> m_mvtxgy = {};
0219   std::vector<float> m_mvtxgz = {};
0220   std::vector<float> m_mvtxgr = {};
0221   float m_projgx = std::numeric_limits<float>::quiet_NaN();
0222   float m_projgy = std::numeric_limits<float>::quiet_NaN();
0223   float m_projgz = std::numeric_limits<float>::quiet_NaN();
0224   float m_projgr = std::numeric_limits<float>::quiet_NaN();
0225   float m_projlx = std::numeric_limits<float>::quiet_NaN();
0226   float m_projlz = std::numeric_limits<float>::quiet_NaN();
0227   float m_clusgx = std::numeric_limits<float>::quiet_NaN();
0228   float m_clusgy = std::numeric_limits<float>::quiet_NaN();
0229   float m_clusgz = std::numeric_limits<float>::quiet_NaN();
0230   float m_clusgr = std::numeric_limits<float>::quiet_NaN();
0231   float m_cluslx = std::numeric_limits<float>::quiet_NaN();
0232   float m_cluslz = std::numeric_limits<float>::quiet_NaN();
0233 
0234   ActsGeometry *m_tGeometry = nullptr;
0235   TrackSeedContainer *m_seedContainer = nullptr;
0236   TrkrClusterContainer *m_clusterMap = nullptr;
0237   PHG4CylinderGeomContainer *m_geomContainerIntt = nullptr;
0238 
0239   int m_lowStrobeIndex = 0;
0240   int m_highStrobeIndex = 1;
0241   /// Configuration classes for Acts seeding
0242   Acts::SeedFinderConfig<SpacePoint> m_seedFinderCfg;
0243   Acts::CylindricalSpacePointGridConfig m_gridCfg;
0244   Acts::CylindricalSpacePointGridOptions m_gridOptions;
0245   Acts::SeedFinderOptions m_seedFinderOptions;
0246 
0247   /// boolean whether or not to include the intt in the acts search windows
0248   bool m_searchInIntt = false;
0249 
0250   /// Configurable parameters
0251   /// seed pt has to be in MeV
0252   float m_minSeedPt = 100 * Acts::UnitConstants::MeV;
0253   float m_uncfactor = 3.18;
0254 
0255   /// How many seeds a given hit can be the middle hit of the seed
0256   /// MVTX can only have the middle layer be the middle hit
0257   int m_maxSeedsPerSpM = 1;
0258 
0259   /// Limiting location of measurements (e.g. detector constraints)
0260   float m_rMax = 200. * Acts::UnitConstants::mm;
0261   float m_rMin = 15. * Acts::UnitConstants::mm;
0262   float m_zMax = 500. * Acts::UnitConstants::mm;
0263   float m_zMin = -500. * Acts::UnitConstants::mm;
0264 
0265   /// misalignment parameters
0266   float m_helixcut = 1;
0267   float m_tolerance = 1.1 * Acts::UnitConstants::mm;
0268   float m_ralign = 0;
0269   float m_zalign = 0;
0270   float m_maxPtScattering = 10;
0271   float m_sigmaScattering = 5.;
0272   float m_sigmaError = 5;
0273 
0274   /// Value tuned to provide as large of phi bins as possible.
0275   /// Increases the secondary finding efficiency
0276   float m_gridFactor = 2.3809;
0277 
0278   /// max distance between two measurements in one seed
0279   float m_deltaRMax = 15 * Acts::UnitConstants::mm;
0280   float m_deltaRMin = 1. * Acts::UnitConstants::mm;
0281   /// Cot of maximum theta angle
0282   float m_cotThetaMax = 2.9;
0283 
0284   /// Maximum impact parameter allowed in mm
0285   float m_impactMax = 20 * Acts::UnitConstants::mm;
0286 
0287   /// Only used in seeding with specified z bin edges, which
0288   /// is more configuration than we need
0289   int m_numPhiNeighbors = 1;
0290 
0291   /// B field value in z direction
0292   /// bfield for space point grid neds to be in kiloTesla
0293   float m_bField = 1.4 * Acts::UnitConstants::T;
0294   std::vector<std::pair<int, int>> zBinNeighborsTop;
0295   std::vector<std::pair<int, int>> zBinNeighborsBottom;
0296   int nphineighbors = 1;
0297  std::unique_ptr<const Acts::GridBinFinder<2ul>> m_bottomBinFinder;
0298  std::unique_ptr<const Acts::GridBinFinder<2ul>> m_topBinFinder;
0299 
0300  int m_event = 0;
0301 
0302  /// Maximum allowed transverse PCA for seed, cm
0303  double m_maxSeedPCA = 2.;
0304 
0305 
0306   /// Search window for phi to match intt clusters in cm
0307   double m_inttrPhiSearchWin = 0.1;
0308   float m_inttzSearchWin = 0.8;  // default to a half strip width
0309   double m_mvtxrPhiSearchWin = 0.2;
0310   float m_mvtxzSearchWin = 0.5;
0311   /// Whether or not to use truth clusters in hit lookup
0312   bool m_useTruthClusters = false;
0313 
0314  bool m_cleanSeeds = false;
0315 
0316  int m_nBadUpdates = 0;
0317  int m_nBadInitialFits = 0;
0318  TrkrClusterIterationMapv1 *_iteration_map = nullptr;
0319  int m_nIteration = 0;
0320  std::string _track_map_name = "SiliconTrackSeedContainer";
0321  ClusterErrorPara _ClusErrPara;
0322 
0323  bool m_seedAnalysis = false;
0324  TFile *m_file = nullptr;
0325  TH2 *h_nInttProj = nullptr;
0326  TH1 *h_nMvtxHits = nullptr;
0327  TH1 *h_nInttHits = nullptr;
0328  TH1 *h_nMatchedClusters = nullptr;
0329  TH2 *h_nHits = nullptr;
0330  TH1 *h_nSeeds = nullptr;
0331  TH1 *h_nActsSeeds = nullptr;
0332  TH1 *h_nTotSeeds = nullptr;
0333  TH1 *h_nInputMeas = nullptr;
0334  TH1 *h_nInputMvtxMeas = nullptr;
0335  TH1 *h_nInputInttMeas = nullptr;
0336  TH2 *h_hits = nullptr;
0337  TH2 *h_zhits = nullptr;
0338  TH2 *h_projHits = nullptr;
0339  TH2 *h_zprojHits = nullptr;
0340  TH2 *h_resids = nullptr;
0341 };
0342 
0343 #endif