macros/common/Trkr_Reco.C

0001 #ifndef MACRO_TRKRRECO_C
0002 #define MACRO_TRKRRECO_C
0003
0004 #include <GlobalVariables.C>
0005
0006 #include <G4_TrkrVariables.C>
0007
0008 #include <trackingdiagnostics/TrackSeedTrackMapConverter.h>
0009
0010 #include <trackreco/MakeActsGeometry.h>
0011 #include <trackreco/PHActsSiliconSeeding.h>
0012 #include <trackreco/PHActsTrackProjection.h>
0013 #include <trackreco/PHActsTrkFitter.h>
0014 #include <trackreco/PHActsVertexPropagator.h>
0015 #include <trackreco/PHCASeeding.h>
0016 #include <trackreco/PHGenFitTrkFitter.h>
0017 #include <trackreco/PHMicromegasTpcTrackMatching.h>
0018 #include <trackreco/PHSiliconHelicalPropagator.h>
0019 #include <trackreco/PHSiliconSeedMerger.h>
0020 #include <trackreco/PHSiliconTpcTrackMatching.h>
0021 #include <trackreco/PHSimpleKFProp.h>
0022 #include <trackreco/PHSimpleVertexFinder.h>
0023 #include <trackreco/PHTpcDeltaZCorrection.h>
0024 #include <trackreco/PHTrackCleaner.h>
0025 #include <trackreco/PHTrackSeeding.h>
0026 #include <trackreco/PrelimDistortionCorrection.h>
0027 #include <trackreco/SecondaryVertexFinder.h>
0028 #include <trackreco/TrackingIterationCounter.h>
0029
0030 #include <tpc/LaserEventRejecter.h>
0031 #include <tpc/TpcLoadDistortionCorrection.h>
0032
0033 #include <tpccalib/PHTpcResiduals.h>
0034
0035 #include <trackermillepedealignment/HelicalFitter.h>
0036 #include <trackermillepedealignment/MakeMilleFiles.h>
0037
0038 #include <fun4all/Fun4AllServer.h>
0039
0040 #include <trackingdiagnostics/TrackContainerCombiner.h>
0041
0042 #include <string>
0043
0044 R__LOAD_LIBRARY(libTrackingDiagnostics.so)
0045 R__LOAD_LIBRARY(libtrack_reco.so)
0046 R__LOAD_LIBRARY(libtpccalib.so)
0047 R__LOAD_LIBRARY(libtpc.so)
0048 R__LOAD_LIBRARY(libtrackeralign.so)
0049
0050 void convert_seeds()
0051 {
0052   auto *se = Fun4AllServer::instance();
0053   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0054
0055   auto *converter = new TrackSeedTrackMapConverter;
0056   // Default set to full SvtxTrackSeeds. Can be set to
0057   // SiliconTrackSeedContainer or TpcTrackSeedContainer
0058   converter->setTrackSeedName("SvtxTrackSeedContainer");
0059   converter->setFieldMap(G4MAGNET::magfield_tracking);
0060   converter->Verbosity(verbosity);
0061   se->registerSubsystem(converter);
0062 }
0063 void Tracking_Reco_Vertex_run2pp()
0064 {
0065   auto *se = Fun4AllServer::instance();
0066   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0067
0068   auto *finder = new PHSimpleVertexFinder;
0069   finder->Verbosity(verbosity);
0070
0071   // new cuts
0072   finder->setDcaCut(0.05);
0073   finder->setTrackPtCut(0.1);
0074   finder->setBeamLineCut(1);
0075   finder->setTrackQualityCut(300);
0076   finder->setNmvtxRequired(3);
0077   finder->setOutlierPairCut(0.10);
0078
0079   se->registerSubsystem(finder);
0080
0081   // Propagate track positions to the vertex position
0082   auto *vtxProp = new PHActsVertexPropagator;
0083   vtxProp->Verbosity(verbosity);
0084   vtxProp->fieldMap(G4MAGNET::magfield_tracking);
0085   se->registerSubsystem(vtxProp);
0086
0087 }
0088 void Tracking_Reco_TrackFit_run2pp(const std::string &outfile = "run2pptrackfit.root")
0089 {
0090   auto *se = Fun4AllServer::instance();
0091   //  int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0092   auto *deltazcorr = new PHTpcDeltaZCorrection;
0093   deltazcorr->Verbosity(0);
0094   se->registerSubsystem(deltazcorr);
0095
0096   // perform final track fit with ACTS
0097   auto *actsFit = new PHActsTrkFitter;
0098   actsFit->Verbosity(0);
0099   actsFit->commissioning(G4TRACKING::use_alignment);
0100   // in calibration mode, fit only Silicons and Micromegas hits
0101   actsFit->fitSiliconMMs(G4TRACKING::SC_CALIBMODE);
0102   actsFit->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0103   actsFit->set_pp_mode(TRACKING::pp_mode);
0104   actsFit->set_use_clustermover(true);  // default is true for now
0105   actsFit->useActsEvaluator(false);
0106   actsFit->useOutlierFinder(false);
0107   actsFit->setFieldMap(G4MAGNET::magfield_tracking);
0108   se->registerSubsystem(actsFit);
0109
0110   auto *cleaner = new PHTrackCleaner();
0111   cleaner->Verbosity(0);
0112   cleaner->set_pp_mode(TRACKING::pp_mode);
0113   se->registerSubsystem(cleaner);
0114
0115   if (G4TRACKING::SC_CALIBMODE)
0116   {
0117     /*
0118      * in calibration mode, calculate residuals between TPC and fitted tracks,
0119      * store in dedicated structure for distortion correction
0120      */
0121     auto *residuals = new PHTpcResiduals;
0122     const TString tpc_residoutfile = outfile + "_PhTpcResiduals.root";
0123     residuals->setOutputfile(tpc_residoutfile.Data());
0124     residuals->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0125
0126     // matches Tony's analysis
0127     residuals->setMinPt(0.2);
0128
0129     // reconstructed distortion grid size (phi, r, z)
0130     residuals->setGridDimensions(36, 48, 80);
0131     se->registerSubsystem(residuals);
0132   }
0133 }
0134 void Tracking_Reco_TpcSeed_run2pp()
0135 {
0136   auto *se = Fun4AllServer::instance();
0137   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0138
0139   /*
0140    * Tpc Seeding
0141    */
0142   auto *seeder = new PHCASeeding("PHCASeeding");
0143   double fieldstrength = std::numeric_limits<double>::quiet_NaN();  // set by isConstantField if constant
0144   bool ConstField = isConstantField(G4MAGNET::magfield_tracking, fieldstrength);
0145   if (ConstField)
0146   {
0147     seeder->useConstBField(true);
0148     seeder->constBField(fieldstrength);
0149   }
0150   else
0151   {
0152     seeder->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0153     seeder->useConstBField(false);
0154     seeder->magFieldFile(G4MAGNET::magfield_tracking);  // to get charge sign right
0155   }
0156   seeder->Verbosity(verbosity);
0157   seeder->SetLayerRange(7, 55);
0158   seeder->SetSearchWindow(2., 0.05);           // z-width and phi-width, default in macro at 1.5 and 0.05
0159   seeder->SetClusAdd_delta_window(3.0, 0.06);  //  (0.5, 0.005) are default; sdzdr_cutoff, d2/dr2(phi)_cutoff
0160   // seeder->SetNClustersPerSeedRange(4,60); // default is 6, 6
0161   seeder->SetMinHitsPerCluster(0);
0162   seeder->SetMinClustersPerTrack(3);
0163   seeder->useFixedClusterError(true);
0164   seeder->set_pp_mode(true);
0165   seeder->reject_zsize1_clusters(true);
0166   se->registerSubsystem(seeder);
0167
0168   // expand stubs in the TPC using simple kalman filter
0169   auto *cprop = new PHSimpleKFProp("PHSimpleKFProp");
0170   cprop->set_field_dir(G4MAGNET::magfield_rescale);
0171   if (ConstField)
0172   {
0173     cprop->useConstBField(true);
0174     cprop->setConstBField(fieldstrength);
0175   }
0176   else
0177   {
0178     cprop->magFieldFile(G4MAGNET::magfield_tracking);
0179     cprop->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0180   }
0181   cprop->useFixedClusterError(true);
0182   cprop->set_max_window(5.);
0183   cprop->Verbosity(verbosity);
0184   cprop->set_pp_mode(true);
0185   cprop->set_max_seeds(5000);
0186   se->registerSubsystem(cprop);
0187
0188   // Always apply preliminary distortion corrections to TPC clusters before silicon matching
0189   // and refit the trackseeds. Replace KFProp fits with the new fit parameters in the TPC seeds.
0190   auto *prelim_distcorr = new PrelimDistortionCorrection;
0191   prelim_distcorr->set_pp_mode(true);
0192   prelim_distcorr->Verbosity(verbosity);
0193   se->registerSubsystem(prelim_distcorr);
0194 }
0195
0196 void Tracking_Reco_SiliconSeed_run2pp()
0197 {
0198   auto *se = Fun4AllServer::instance();
0199   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0200   /*
0201    * Silicon Seeding
0202    */
0203
0204   auto *silicon_Seeding = new PHActsSiliconSeeding;
0205   silicon_Seeding->Verbosity(verbosity);
0206   silicon_Seeding->setStrobeRange(-5, 5);
0207   silicon_Seeding->isStreaming();
0208   // these get us to about 83% INTT > 1
0209   silicon_Seeding->setinttRPhiSearchWindow(0.2);
0210   se->registerSubsystem(silicon_Seeding);
0211
0212   auto *merger = new PHSiliconSeedMerger;
0213   merger->Verbosity(verbosity);
0214   se->registerSubsystem(merger);
0215 }
0216 void Tracking_Reco_TrackSeed_run2pp()
0217 {
0218   Tracking_Reco_SiliconSeed_run2pp();
0219   Tracking_Reco_TpcSeed_run2pp();
0220 }
0221
0222 void Tracking_Reco_TrackMatching_run2pp()
0223 {
0224   auto *se = Fun4AllServer::instance();
0225   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0226   /*
0227    * Track Matching between silicon and TPC
0228    */
0229   // The normal silicon association methods
0230   // Match the TPC track stubs from the CA seeder to silicon track stubs from PHSiliconTruthTrackSeeding
0231   auto *silicon_match = new PHSiliconTpcTrackMatching;
0232   silicon_match->Verbosity(verbosity);
0233   silicon_match->set_pp_mode(TRACKING::pp_mode);
0234   if (G4TPC::ENABLE_AVERAGE_CORRECTIONS)
0235   {
0236     // for general tracking
0237     // Eta/Phi window is determined by 3 sigma window
0238     // X/Y/Z window is determined by 4 sigma window
0239     silicon_match->window_deta.set_posQoverpT_maxabs({-0.014, 0.0331, 0.48});
0240     silicon_match->window_deta.set_negQoverpT_maxabs({-0.006, 0.0235, 0.52});
0241     silicon_match->set_deltaeta_min(0.03);
0242     silicon_match->window_dphi.set_QoverpT_range({-0.15, 0, 0}, {0.15, 0, 0});
0243     silicon_match->window_dx.set_QoverpT_maxabs({3.0, 0, 0});
0244     silicon_match->window_dy.set_QoverpT_maxabs({3.0, 0, 0});
0245     silicon_match->window_dz.set_posQoverpT_maxabs({1.138, 0.3919, 0.84});
0246     silicon_match->window_dz.set_negQoverpT_maxabs({0.719, 0.6485, 0.65});
0247     silicon_match->set_crossing_deltaz_max(30);
0248     silicon_match->set_crossing_deltaz_min(2);
0249
0250     // for distortion correction using SI-TPOT fit and track pT>0.5
0251     if (G4TRACKING::SC_CALIBMODE)
0252     {
0253       silicon_match->window_deta.set_posQoverpT_maxabs({0.016, 0.0060, 1.13});
0254       silicon_match->window_deta.set_negQoverpT_maxabs({0.022, 0.0022, 1.44});
0255       silicon_match->set_deltaeta_min(0.03);
0256       silicon_match->window_dphi.set_QoverpT_range({-0.15, 0, 0}, {0.09, 0, 0});
0257       silicon_match->window_dx.set_QoverpT_maxabs({2.0, 0, 0});
0258       silicon_match->window_dy.set_QoverpT_maxabs({1.5, 0, 0});
0259       silicon_match->window_dz.set_posQoverpT_maxabs({1.213, 0.0211, 2.09});
0260       silicon_match->window_dz.set_negQoverpT_maxabs({1.307, 0.0001, 4.52});
0261       silicon_match->set_crossing_deltaz_min(1.2);
0262     }
0263   }
0264   se->registerSubsystem(silicon_match);
0265
0266   // Match TPC track stubs from CA seeder to clusters in the micromegas layers
0267   auto *mm_match = new PHMicromegasTpcTrackMatching;
0268   mm_match->Verbosity(verbosity);
0269   mm_match->set_pp_mode(TRACKING::pp_mode);
0270   mm_match->set_rphi_search_window_lyr1(3.);
0271   mm_match->set_rphi_search_window_lyr2(15.0);
0272   mm_match->set_z_search_window_lyr1(30.0);
0273   mm_match->set_z_search_window_lyr2(3.);
0274
0275   mm_match->set_min_tpc_layer(38);             // layer in TPC to start projection fit
0276   mm_match->set_test_windows_printout(false);  // used for tuning search windows only
0277   se->registerSubsystem(mm_match);
0278 }
0279 void Tracking_Reco_TrackSeed_ZeroField()
0280 {
0281   // set up verbosity
0282   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0283
0284   // get fun4all server instance
0285   auto *se = Fun4AllServer::instance();
0286
0287   auto *silicon_Seeding = new PHActsSiliconSeeding;
0288   silicon_Seeding->Verbosity(verbosity);
0289   silicon_Seeding->setStrobeRange(-5, 5);
0290   silicon_Seeding->isStreaming();
0291   // these get us to about 83% INTT > 1
0292   silicon_Seeding->setinttRPhiSearchWindow(0.2);
0293   se->registerSubsystem(silicon_Seeding);
0294
0295   auto *merger = new PHSiliconSeedMerger;
0296   merger->Verbosity(verbosity);
0297   se->registerSubsystem(merger);
0298
0299   auto *seeder = new PHCASeeding("PHCASeeding");
0300   double fieldstrength = std::numeric_limits<double>::quiet_NaN();  // set by isConstantField if constant
0301   bool ConstField = isConstantField(G4MAGNET::magfield_tracking, fieldstrength);
0302   if (ConstField)
0303   {
0304     seeder->useConstBField(true);
0305     seeder->constBField(fieldstrength);
0306   }
0307   else
0308   {
0309     seeder->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0310     seeder->useConstBField(false);
0311     seeder->magFieldFile(G4MAGNET::magfield_tracking);  // to get charge sign right
0312   }
0313   seeder->Verbosity(verbosity);
0314   seeder->SetLayerRange(7, 55);
0315   seeder->SetSearchWindow(1.5, 0.05);  // (z width, phi width)
0316   seeder->SetMinHitsPerCluster(0);
0317   seeder->SetMinClustersPerTrack(3);
0318   seeder->useFixedClusterError(true);
0319
0320   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0321   {
0322     seeder->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0323     seeder->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0324     seeder->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0325     seeder->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0326     seeder->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0327   }
0328   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0329   {
0330     seeder->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0331     seeder->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0332     seeder->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0333     seeder->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0334     seeder->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0335   }
0336   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0337   {
0338     seeder->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0339     seeder->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0340     seeder->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0341     seeder->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0342     seeder->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0343   }
0344   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0345   {
0346     seeder->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0347     seeder->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0348     seeder->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0349     seeder->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0350     seeder->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0351   }
0352   else
0353   {
0354   }
0355
0356   seeder->set_pp_mode(true);
0357   se->registerSubsystem(seeder);
0358
0359   // expand stubs in the TPC using simple kalman filter
0360   auto *cprop = new PHSimpleKFProp("PHSimpleKFProp");
0361   cprop->set_field_dir(G4MAGNET::magfield_rescale);
0362   if (ConstField)
0363   {
0364     cprop->useConstBField(true);
0365     cprop->setConstBField(fieldstrength);
0366   }
0367   else
0368   {
0369     cprop->magFieldFile(G4MAGNET::magfield_tracking);
0370     cprop->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0371   }
0372   cprop->useFixedClusterError(true);
0373   cprop->set_max_window(5.);
0374   cprop->Verbosity(verbosity);
0375
0376   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0377   {
0378     cprop->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0379     cprop->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0380     cprop->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0381     cprop->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0382     cprop->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0383   }
0384   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0385   {
0386     cprop->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0387     cprop->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0388     cprop->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0389     cprop->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0390     cprop->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0391   }
0392   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0393   {
0394     cprop->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0395     cprop->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0396     cprop->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0397     cprop->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0398     cprop->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0399   }
0400   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0401   {
0402     cprop->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0403     cprop->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0404     cprop->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0405     cprop->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0406     cprop->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0407   }
0408   else
0409   {
0410   }
0411
0412   cprop->set_pp_mode(true);
0413   se->registerSubsystem(cprop);
0414 }
0415 void Tracking_Reco_TrackSeed()
0416 {
0417   // set up verbosity
0418   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0419
0420   // get fun4all server instance
0421   auto *se = Fun4AllServer::instance();
0422
0423   // Assemble silicon clusters into track stubs
0424
0425   auto *silicon_Seeding = new PHActsSiliconSeeding;
0426   silicon_Seeding->Verbosity(verbosity);
0427   silicon_Seeding->isStreaming();
0428   // modify strobe range
0429   silicon_Seeding->setStrobeRange(-1, 2);
0430
0431   se->registerSubsystem(silicon_Seeding);
0432
0433   auto *merger = new PHSiliconSeedMerger;
0434   merger->Verbosity(verbosity);
0435   se->registerSubsystem(merger);
0436
0437   auto *seeder = new PHCASeeding("PHCASeeding");
0438   double fieldstrength = std::numeric_limits<double>::quiet_NaN();  // set by isConstantField if constant
0439   bool ConstField = isConstantField(G4MAGNET::magfield_tracking, fieldstrength);
0440   if (ConstField)
0441   {
0442     seeder->useConstBField(true);
0443     seeder->constBField(fieldstrength);
0444   }
0445   else
0446   {
0447     seeder->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0448     seeder->useConstBField(false);
0449     seeder->magFieldFile(G4MAGNET::magfield_tracking);  // to get charge sign right
0450   }
0451   seeder->Verbosity(verbosity);
0452   seeder->SetLayerRange(7, 55);
0453   seeder->SetSearchWindow(1.5, 0.05);  // (z width, phi width)
0454   seeder->SetMinHitsPerCluster(0);
0455   seeder->SetMinClustersPerTrack(3);
0456   seeder->useFixedClusterError(true);
0457
0458   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0459   {
0460     seeder->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0461     seeder->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0462     seeder->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0463     seeder->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0464     seeder->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0465   }
0466   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0467   {
0468     seeder->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0469     seeder->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0470     seeder->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0471     seeder->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0472     seeder->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0473   }
0474   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0475   {
0476     seeder->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0477     seeder->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0478     seeder->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0479     seeder->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0480     seeder->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0481   }
0482   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0483   {
0484     seeder->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0485     seeder->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0486     seeder->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0487     seeder->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0488     seeder->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0489   }
0490   else
0491   {
0492   }
0493
0494   seeder->set_pp_mode(TRACKING::pp_mode);
0495   se->registerSubsystem(seeder);
0496
0497   // expand stubs in the TPC using simple kalman filter
0498   auto *cprop = new PHSimpleKFProp("PHSimpleKFProp");
0499   cprop->set_field_dir(G4MAGNET::magfield_rescale);
0500   if (ConstField)
0501   {
0502     cprop->useConstBField(true);
0503     cprop->setConstBField(fieldstrength);
0504   }
0505   else
0506   {
0507     cprop->magFieldFile(G4MAGNET::magfield_tracking);
0508     cprop->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0509   }
0510   cprop->useFixedClusterError(true);
0511   cprop->set_max_window(5.);
0512   cprop->Verbosity(verbosity);
0513
0514   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0515   {
0516     cprop->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0517     cprop->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0518     cprop->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0519     cprop->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0520     cprop->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0521   }
0522   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0523   {
0524     cprop->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0525     cprop->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0526     cprop->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0527     cprop->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0528     cprop->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0529   }
0530   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0531   {
0532     cprop->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0533     cprop->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0534     cprop->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0535     cprop->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0536     cprop->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0537   }
0538   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0539   {
0540     cprop->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0541     cprop->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0542     cprop->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0543     cprop->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0544     cprop->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0545   }
0546   else
0547   {
0548   }
0549
0550   cprop->set_pp_mode(TRACKING::pp_mode);
0551   se->registerSubsystem(cprop);
0552
0553   if (TRACKING::pp_mode)
0554   {
0555     // for pp mode, apply preliminary distortion corrections to TPC clusters before crossing is known
0556     // and refit the trackseeds. Replace KFProp fits with the new fit parameters in the TPC seeds.
0557     auto *prelim_distcorr = new PrelimDistortionCorrection;
0558     prelim_distcorr->set_pp_mode(TRACKING::pp_mode);
0559     prelim_distcorr->Verbosity(verbosity);
0560
0561     if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0562     {
0563       prelim_distcorr->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0564       prelim_distcorr->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0565       prelim_distcorr->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0566       prelim_distcorr->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0567       prelim_distcorr->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0568     }
0569     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0570     {
0571       prelim_distcorr->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0572       prelim_distcorr->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0573       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0574       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0575       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0576     }
0577     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0578     {
0579       prelim_distcorr->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0580       prelim_distcorr->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0581       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0582       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0583       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0584     }
0585     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0586     {
0587       prelim_distcorr->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0588       prelim_distcorr->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0589       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0590       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0591       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0592     }
0593     else
0594     {
0595     }
0596
0597     se->registerSubsystem(prelim_distcorr);
0598   }
0599
0600   std::cout << "Tracking_Reco_TrackSeed - Using stub matching for Si matching " << std::endl;
0601   // The normal silicon association methods
0602   // Match the TPC track stubs from the CA seeder to silicon track stubs from PHSiliconTruthTrackSeeding
0603   auto *silicon_match = new PHSiliconTpcTrackMatching;
0604   silicon_match->Verbosity(verbosity);
0605   silicon_match->set_pp_mode(TRACKING::pp_mode);
0606   std::cout << "PHSiliconTpcTrackMatching pp_mode set to " << TRACKING::pp_mode << std::endl;
0607   if (G4TRACKING::SC_CALIBMODE)
0608   {
0609     // search windows for initial matching with distortions
0610     // tuned values are 0.04 and 0.008 in distorted events
0611     silicon_match->set_phi_search_window(0.04);
0612     silicon_match->set_eta_search_window(0.008);
0613   }
0614   else
0615   {
0616     // after distortion corrections and rerunning clustering, default tuned values are 0.02 and 0.004 in low occupancy events
0617     silicon_match->set_phi_search_window(0.03);
0618     silicon_match->set_eta_search_window(0.005);
0619   }
0620   silicon_match->set_test_windows_printout(false);  // used for tuning search windows
0621   se->registerSubsystem(silicon_match);
0622
0623   // Associate Micromegas clusters with the tracks
0624   if (Enable::MICROMEGAS)
0625   {
0626     std::cout << "Tracking_Reco_TrackSeed - Using Micromegas matching " << std::endl;
0627
0628     // Match TPC track stubs from CA seeder to clusters in the micromegas layers
0629     auto *mm_match = new PHMicromegasTpcTrackMatching;
0630     mm_match->Verbosity(verbosity);
0631     mm_match->set_pp_mode(TRACKING::pp_mode);
0632     mm_match->set_rphi_search_window_lyr1(0.2);
0633     mm_match->set_rphi_search_window_lyr2(13.0);
0634     mm_match->set_z_search_window_lyr1(26.0);
0635     mm_match->set_z_search_window_lyr2(0.2);
0636
0637     mm_match->set_min_tpc_layer(38);             // layer in TPC to start projection fit
0638     mm_match->set_test_windows_printout(false);  // used for tuning search windows only
0639     se->registerSubsystem(mm_match);
0640   }
0641 }
0642
0643 void Tracking_Reco_TrackSeed_pass1()
0644 {
0645   Fun4AllServer *se = Fun4AllServer::instance();
0646   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0647
0648   TrackingIterationCounter *counter = new TrackingIterationCounter("TrkrIter1");
0649   /// Clusters already used are in the 0th iteration
0650   counter->iteration(0);
0651   se->registerSubsystem(counter);
0652
0653   PHActsSiliconSeeding *silseed = new PHActsSiliconSeeding("PHActsSiliconSeedingIt1");
0654   silseed->Verbosity(verbosity);
0655   silseed->searchInIntt();
0656   silseed->iteration(1);
0657   silseed->set_track_map_name("SiliconTrackSeedContainerIt1");
0658   se->registerSubsystem(silseed);
0659
0660   PHSiliconSeedMerger *merger = new PHSiliconSeedMerger("SiliconSeedMargerIt1");
0661   merger->Verbosity(verbosity);
0662   merger->clusterOverlap(2);
0663   merger->searchIntt();
0664   merger->trackMapName("SiliconTrackSeedContainerIt1");
0665   se->registerSubsystem(merger);
0666
0667   TrackContainerCombiner *combiner = new TrackContainerCombiner;
0668   combiner->newContainerName("SiliconTrackSeedContainer");
0669   combiner->oldContainerName("SiliconTrackSeedContainerIt1");
0670   combiner->Verbosity(verbosity);
0671   se->registerSubsystem(combiner);
0672 }
0673
0674 void vertexing()
0675 {
0676   Fun4AllServer *se = Fun4AllServer::instance();
0677   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0678
0679   auto *vtxfinder = new PHSimpleVertexFinder;
0680   vtxfinder->Verbosity(verbosity);
0681   se->registerSubsystem(vtxfinder);
0682 }
0683
0684 void Tracking_Reco_TrackFit()
0685 {
0686   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0687   auto *se = Fun4AllServer::instance();
0688
0689   // correct clusters for particle propagation in TPC
0690   auto *deltazcorr = new PHTpcDeltaZCorrection;
0691   deltazcorr->Verbosity(verbosity);
0692   se->registerSubsystem(deltazcorr);
0693
0694   if (G4TRACKING::use_genfit_track_fitter)
0695   {
0696     // perform final track fit with GENFIT
0697     auto *genfitFit = new PHGenFitTrkFitter;
0698     genfitFit->Verbosity(verbosity);
0699     genfitFit->set_fit_silicon_mms(G4TRACKING::SC_CALIBMODE);
0700     se->registerSubsystem(genfitFit);
0701
0702     if (G4TRACKING::SC_CALIBMODE)
0703     {
0704
0705       /*
0706        * in calibration mode, calculate residuals between TPC and fitted tracks,
0707        * store in dedicated structure for distortion correction
0708        */
0709       auto *residuals = new PHTpcResiduals();
0710       residuals->setTrackMapName("SvtxTrackMap");
0711       residuals->setOutputfile(G4TRACKING::SC_ROOTOUTPUT_FILENAME);
0712       residuals->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0713       // reconstructed distortion grid size (phi, r, z)
0714       residuals->setGridDimensions(36, 48, 80);
0715       residuals->Verbosity(verbosity);
0716       se->registerSubsystem(residuals);
0717     }
0718
0719   }
0720   else
0721   {
0722     // perform final track fit with ACTS
0723     auto *actsFit = new PHActsTrkFitter;
0724     actsFit->Verbosity(verbosity);
0725     actsFit->commissioning(G4TRACKING::use_alignment);
0726     // in calibration mode, fit only Silicons and Micromegas hits
0727     actsFit->fitSiliconMMs(G4TRACKING::SC_CALIBMODE);
0728     actsFit->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0729     actsFit->set_pp_mode(TRACKING::pp_mode);
0730     actsFit->set_use_clustermover(true);  // default is true for now
0731     actsFit->useActsEvaluator(false);
0732     actsFit->useOutlierFinder(false);
0733     actsFit->setFieldMap(G4MAGNET::magfield_tracking);
0734     se->registerSubsystem(actsFit);
0735
0736     if (G4TRACKING::SC_CALIBMODE)
0737     {
0738       /*
0739        * in calibration mode, calculate residuals between TPC and fitted tracks,
0740        * store in dedicated structure for distortion correction
0741        */
0742       auto *residuals = new PHTpcResiduals();
0743       residuals->setOutputfile(G4TRACKING::SC_ROOTOUTPUT_FILENAME);
0744       residuals->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0745       // reconstructed distortion grid size (phi, r, z)
0746       residuals->setGridDimensions(36, 48, 80);
0747       residuals->Verbosity(verbosity);
0748       se->registerSubsystem(residuals);
0749     }
0750   }
0751
0752   if (!G4TRACKING::SC_CALIBMODE)
0753   {
0754     /*
0755      * in full tracking mode, run track cleaner, vertex finder,
0756      * propagete tracks to vertex
0757      * propagate tracks to EMCAL
0758      */
0759
0760     if (!G4TRACKING::use_full_truth_track_seeding)
0761     {
0762       // Choose the best silicon matched track for each TPC track seed
0763       /* this breaks in truth_track seeding mode because there is no TpcSeed */
0764       auto *cleaner = new PHTrackCleaner;
0765       cleaner->Verbosity(verbosity);
0766       // cleaner->set_quality_cut(30.0);
0767       se->registerSubsystem(cleaner);
0768     }
0769
0770     vertexing();
0771
0772     if (!G4TRACKING::use_genfit_track_fitter)
0773     {
0774       // Propagate track positions to the vertex position
0775       auto *vtxProp = new PHActsVertexPropagator;
0776       vtxProp->Verbosity(verbosity);
0777       vtxProp->fieldMap(G4MAGNET::magfield_tracking);
0778       se->registerSubsystem(vtxProp);
0779
0780       // project tracks to EMCAL
0781       auto *projection = new PHActsTrackProjection;
0782       projection->Verbosity(verbosity);
0783       double fieldstrength = std::numeric_limits<double>::quiet_NaN();
0784       if (isConstantField(G4MAGNET::magfield_tracking, fieldstrength))
0785       {
0786         projection->setConstFieldVal(fieldstrength);
0787       }
0788       se->registerSubsystem(projection);
0789     }
0790   }
0791 }
0792
0793 void Tracking_Reco_CommissioningTrackSeed()
0794 {
0795   // set up verbosity
0796   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
0797
0798   // get fun4all server instance
0799   auto *se = Fun4AllServer::instance();
0800
0801   auto *silicon_Seeding = new PHActsSiliconSeeding;
0802   silicon_Seeding->Verbosity(verbosity);
0803   silicon_Seeding->sigmaScattering(50.);
0804   silicon_Seeding->setinttRPhiSearchWindow(2.);
0805   silicon_Seeding->helixcut(0.01);
0806   se->registerSubsystem(silicon_Seeding);
0807
0808   auto *merger = new PHSiliconSeedMerger;
0809   merger->Verbosity(verbosity);
0810   se->registerSubsystem(merger);
0811
0812   // Assemble TPC clusters into track stubs
0813   auto *seeder = new PHCASeeding("PHCASeeding");
0814   seeder->set_field_dir(G4MAGNET::magfield_rescale);                // to get charge sign right
0815   double fieldstrength = std::numeric_limits<double>::quiet_NaN();  // set by isConstantField if constant
0816   bool ConstField = isConstantField(G4MAGNET::magfield_tracking, fieldstrength);
0817   if (!ConstField)
0818   {
0819     seeder->magFieldFile(G4MAGNET::magfield_tracking);
0820     seeder->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0821   }
0822   seeder->Verbosity(verbosity);
0823   seeder->SetLayerRange(7, 55);
0824   seeder->SetSearchWindow(1.5, 0.05);  // (z width, phi width)
0825   seeder->SetMinHitsPerCluster(0);
0826   seeder->SetMinClustersPerTrack(3);
0827   seeder->useConstBField(false);
0828   seeder->useFixedClusterError(true);
0829
0830   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0831   {
0832     seeder->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0833     seeder->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0834     seeder->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0835     seeder->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0836     seeder->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0837   }
0838   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0839   {
0840     seeder->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0841     seeder->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0842     seeder->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0843     seeder->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0844     seeder->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0845   }
0846   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0847   {
0848     seeder->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0849     seeder->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0850     seeder->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0851     seeder->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0852     seeder->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0853   }
0854   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0855   {
0856     seeder->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0857     seeder->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0858     seeder->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0859     seeder->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0860     seeder->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0861   }
0862   else
0863   {
0864   }
0865
0866   seeder->set_pp_mode(TRACKING::pp_mode);
0867   se->registerSubsystem(seeder);
0868
0869   // expand stubs in the TPC using simple kalman filter
0870   auto *cprop = new PHSimpleKFProp("PHSimpleKFProp");
0871   cprop->set_field_dir(G4MAGNET::magfield_rescale);
0872   if (!ConstField)
0873   {
0874     cprop->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0875   }
0876   cprop->useConstBField(false);
0877   cprop->useFixedClusterError(true);
0878   cprop->set_max_window(5.);
0879   cprop->Verbosity(verbosity);
0880
0881   if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0882   {
0883     cprop->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0884     cprop->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0885     cprop->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0886     cprop->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0887     cprop->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0888   }
0889   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0890   {
0891     cprop->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0892     cprop->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0893     cprop->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0894     cprop->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0895     cprop->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0896   }
0897   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0898   {
0899     cprop->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0900     cprop->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0901     cprop->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0902     cprop->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0903     cprop->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0904   }
0905   else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0906   {
0907     cprop->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0908     cprop->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0909     cprop->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0910     cprop->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0911     cprop->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0912   }
0913   else
0914   {
0915   }
0916
0917   cprop->set_pp_mode(TRACKING::pp_mode);
0918   se->registerSubsystem(cprop);
0919
0920   if (TRACKING::pp_mode)
0921   {
0922     // for pp mode, apply preliminary distortion corrections to TPC clusters before crossing is known
0923     // and refit the trackseeds. Replace KFProp fits with the new fit parameters in the TPC seeds.
0924     auto *prelim_distcorr = new PrelimDistortionCorrection;
0925     prelim_distcorr->set_pp_mode(TRACKING::pp_mode);
0926     prelim_distcorr->Verbosity(verbosity);
0927
0928     if (G4TPC::TPC_GAS_MIXTURE == "NeCF4")
0929     {
0930       prelim_distcorr->setNeonFraction(G4TPC::NeCF4_Ne_frac);
0931       prelim_distcorr->setArgonFraction(G4TPC::NeCF4_Ar_frac);
0932       prelim_distcorr->setCF4Fraction(G4TPC::NeCF4_CF4_frac);
0933       prelim_distcorr->setNitrogenFraction(G4TPC::NeCF4_N2_frac);
0934       prelim_distcorr->setIsobutaneFraction(G4TPC::NeCF4_isobutane_frac);
0935     }
0936     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4")
0937     {
0938       prelim_distcorr->setNeonFraction(G4TPC::ArCF4_Ne_frac);
0939       prelim_distcorr->setArgonFraction(G4TPC::ArCF4_Ar_frac);
0940       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4_CF4_frac);
0941       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4_N2_frac);
0942       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4_isobutane_frac);
0943     }
0944     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4N2")
0945     {
0946       prelim_distcorr->setNeonFraction(G4TPC::ArCF4N2_Ne_frac);
0947       prelim_distcorr->setArgonFraction(G4TPC::ArCF4N2_Ar_frac);
0948       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4N2_CF4_frac);
0949       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4N2_N2_frac);
0950       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4N2_isobutane_frac);
0951     }
0952     else if (G4TPC::TPC_GAS_MIXTURE == "ArCF4Isobutane")
0953     {
0954       prelim_distcorr->setNeonFraction(G4TPC::ArCF4Isobutane_Ne_frac);
0955       prelim_distcorr->setArgonFraction(G4TPC::ArCF4Isobutane_Ar_frac);
0956       prelim_distcorr->setCF4Fraction(G4TPC::ArCF4Isobutane_CF4_frac);
0957       prelim_distcorr->setNitrogenFraction(G4TPC::ArCF4Isobutane_N2_frac);
0958       prelim_distcorr->setIsobutaneFraction(G4TPC::ArCF4Isobutane_isobutane_frac);
0959     }
0960     else
0961     {
0962     }
0963
0964     se->registerSubsystem(prelim_distcorr);
0965   }
0966
0967   // match silicon track seeds to TPC track seeds
0968
0969   // The normal silicon association methods
0970   // Match the TPC track stubs from the CA seeder to silicon track stubs from PHSiliconTruthTrackSeeding
0971   auto *silicon_match = new PHSiliconTpcTrackMatching;
0972   silicon_match->Verbosity(verbosity);
0973   silicon_match->set_pp_mode(TRACKING::pp_mode);
0974
0975   silicon_match->set_phi_search_window(0.2);
0976   silicon_match->set_eta_search_window(0.015);
0977   silicon_match->set_x_search_window(std::numeric_limits<double>::max());
0978   silicon_match->set_y_search_window(std::numeric_limits<double>::max());
0979   silicon_match->set_z_search_window(std::numeric_limits<double>::max());
0980
0981   silicon_match->set_test_windows_printout(false);  // used for tuning search windows
0982   se->registerSubsystem(silicon_match);
0983
0984   // Associate Micromegas clusters with the tracks
0985   if (Enable::MICROMEGAS)
0986   {
0987     // Match TPC track stubs from CA seeder to clusters in the micromegas layers
0988     auto *mm_match = new PHMicromegasTpcTrackMatching;
0989     mm_match->Verbosity(verbosity);
0990
0991     mm_match->set_rphi_search_window_lyr1(0.4);
0992     mm_match->set_rphi_search_window_lyr2(13.0);
0993     mm_match->set_z_search_window_lyr1(26.0);
0994     mm_match->set_z_search_window_lyr2(0.2);
0995
0996     mm_match->set_min_tpc_layer(38);             // layer in TPC to start projection fit
0997     mm_match->set_test_windows_printout(false);  // used for tuning search windows only
0998     se->registerSubsystem(mm_match);
0999   }
1000 }
1001
1002 void alignment(const std::string &datafilename = "mille_output_data_file",
1003                const std::string &steeringfilename = "mille_steer")
1004 {
1005   Fun4AllServer *se = Fun4AllServer::instance();
1006   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
1007
1008   auto *mille = new MakeMilleFiles;
1009   mille->Verbosity(verbosity);
1010   mille->set_datafile_name(datafilename + ".bin");
1011   mille->set_steeringfile_name(steeringfilename + ".txt");
1012   se->registerSubsystem(mille);
1013
1014   auto *helical = new HelicalFitter;
1015   helical->Verbosity(verbosity);
1016   helical->set_datafile_name(datafilename + "_helical.bin");
1017   helical->set_steeringfile_name(steeringfilename + "_helical.txt");
1018   se->registerSubsystem(helical);
1019 }
1020
1021 void Tracking_Reco()
1022 {
1023   /*
1024    * just a wrapper around track seeding and track fitting methods,
1025    * to minimize disruption to existing steering macros
1026    */
1027   if (G4TRACKING::use_alignment)
1028   {
1029     Tracking_Reco_CommissioningTrackSeed();
1030   }
1031   else
1032   {
1033     Tracking_Reco_TrackSeed();
1034   }
1035
1036   if (G4TRACKING::convert_seeds_to_svtxtracks)
1037   {
1038     convert_seeds();
1039     vertexing();
1040   }
1041   else
1042   {
1043     Tracking_Reco_TrackFit();
1044   }
1045
1046   if (G4TRACKING::iterative_seeding)
1047   {
1048     Tracking_Reco_TrackSeed_pass1();
1049
1050     if (G4TRACKING::convert_seeds_to_svtxtracks)
1051     {
1052       convert_seeds();
1053       vertexing();
1054     }
1055   }
1056
1057   if (G4TRACKING::use_alignment)
1058   {
1059     alignment();
1060   }
1061 }
1062
1063 void Filter_Conversion_Electrons(const std::string &ntuple_outfile)
1064 {
1065   Fun4AllServer *se = Fun4AllServer::instance();
1066   int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
1067
1068   SecondaryVertexFinder *secvert = new SecondaryVertexFinder;
1069   secvert->Verbosity(verbosity);
1070   //  secvert->set_write_electrons_node(true);  // writes copy of filtered electron tracks to node tree
1071   //  secvert->set_write_ntuple(false);  // writes ntuple for tuning cuts
1072   secvert->setDecayParticleMass(0.000511);  // for electrons
1073   secvert->setOutfileName(ntuple_outfile);
1074   se->registerSubsystem(secvert);
1075 }
1076
1077 void Reject_Laser_Events()
1078 {
1079   if (G4TPC::REJECT_LASER_EVENTS)
1080   {
1081     auto *se = Fun4AllServer::instance();
1082     //    int verbosity = std::max(Enable::VERBOSITY, Enable::TRACKING_VERBOSITY);
1083
1084     LaserEventRejecter *rejecter = new LaserEventRejecter();
1085     se->registerSubsystem(rejecter);
1086   }
1087 }
1088
1089 #endif