Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-12-18 09:22:22 

0001 /*
0002  * This macro shows running the full event combining + tracking for
0003  * cosmics running.. There are some analysis modules run at the end
0004  * which package  hits, clusters, and clusters on tracks into trees
0005  * for analysis.
0006  */
0007 
0008 // leave the GlobalVariables.C at the beginning, an empty line afterwards
0009 // protects its position against reshuffling by clang-format
0010 #include <GlobalVariables.C>
0011 
0012 #include <G4_ActsGeom.C>
0013 #include <G4_Magnet.C>
0014 #include <QA.C>
0015 #include <Trkr_Clustering.C>
0016 #include <Trkr_Reco.C>
0017 #include <Trkr_RecoInit.C>
0018 #include <Trkr_TpcReadoutInit.C>
0019 
0020 #include <trackreco/PHActsSiliconSeeding.h>
0021 #include <trackreco/PHSiliconSeedMerger.h>
0022 
0023 #include <cdbobjects/CDBTTree.h>
0024 
0025 #include <ffamodules/CDBInterface.h>
0026 #include <ffamodules/FlagHandler.h>
0027 #include <ffamodules/HeadReco.h>
0028 #include <ffamodules/SyncReco.h>
0029 
0030 #include <mvtxrawhitqa/MvtxRawHitQA.h>
0031 
0032 #include <inttrawhitqa/InttRawHitQA.h>
0033 
0034 #include <tpcqa/TpcRawHitQA.h>
0035 
0036 #include <trackingqa/InttClusterQA.h>
0037 #include <trackingqa/MicromegasClusterQA.h>
0038 #include <trackingqa/MvtxClusterQA.h>
0039 #include <trackingqa/SiliconSeedsQA.h>
0040 #include <trackingqa/TpcClusterQA.h>
0041 #include <trackingqa/TpcSeedsQA.h>
0042 #include <trackingqa/TpcSiliconQA.h>
0043 
0044 #include <trackingdiagnostics/TrackResiduals.h>
0045 #include <trackingdiagnostics/TrkrNtuplizer.h>
0046 
0047 #include <fun4allraw/Fun4AllStreamingInputManager.h>
0048 #include <fun4allraw/InputManagerType.h>
0049 #include <fun4allraw/SingleGl1PoolInput.h>
0050 #include <fun4allraw/SingleInttPoolInput.h>
0051 #include <fun4allraw/SingleMicromegasPoolInput.h>
0052 #include <fun4allraw/SingleMvtxPoolInput.h>
0053 #include <fun4allraw/SingleTpcTimeFrameInput.h>
0054 
0055 #include <fun4all/Fun4AllDstInputManager.h>
0056 #include <fun4all/Fun4AllDstOutputManager.h>
0057 #include <fun4all/Fun4AllInputManager.h>
0058 #include <fun4all/Fun4AllOutputManager.h>
0059 #include <fun4all/Fun4AllRunNodeInputManager.h>
0060 #include <fun4all/Fun4AllServer.h>
0061 #include <fun4all/Fun4AllUtils.h>
0062 
0063 #include <phool/recoConsts.h>
0064 
0065 R__LOAD_LIBRARY(libfun4all.so)
0066 R__LOAD_LIBRARY(libffamodules.so)
0067 R__LOAD_LIBRARY(libphool.so)
0068 R__LOAD_LIBRARY(libcdbobjects.so)
0069 R__LOAD_LIBRARY(libmvtx.so)
0070 R__LOAD_LIBRARY(libintt.so)
0071 R__LOAD_LIBRARY(libtpc.so)
0072 R__LOAD_LIBRARY(libmicromegas.so)
0073 R__LOAD_LIBRARY(libTrackingDiagnostics.so)
0074 R__LOAD_LIBRARY(libtrackingqa.so)
0075 R__LOAD_LIBRARY(libtrack_reco.so)
0076 
0077 bool isGood(const std::string &infile);
0078 int getrunnumber(const std::string &listfile);
0079 
0080 void Fun4All_PRDFReconstruction(
0081     const int nEvents = 0,
0082     const int  /*runnumber_unused*/ = 53756, // just kept to preserve the API
0083     const std::string & /*dir*/ = "/.",
0084     const std::string &outfilename = "output_tracks",
0085     const std::string &input_gl1file = "gl1daq.list",
0086     const std::string &input_inttfile00 = "intt0.list",
0087     const std::string &input_inttfile01 = "intt1.list",
0088     const std::string &input_inttfile02 = "intt2.list",
0089     const std::string &input_inttfile03 = "intt3.list",
0090     const std::string &input_inttfile04 = "intt4.list",
0091     const std::string &input_inttfile05 = "intt5.list",
0092     const std::string &input_inttfile06 = "intt6.list",
0093     const std::string &input_inttfile07 = "intt7.list",
0094     const std::string &input_mvtxfile00 = "mvtx0.list",
0095     const std::string &input_mvtxfile01 = "mvtx1.list",
0096     const std::string &input_mvtxfile02 = "mvtx2.list",
0097     const std::string &input_mvtxfile03 = "mvtx3.list",
0098     const std::string &input_mvtxfile04 = "mvtx4.list",
0099     const std::string &input_mvtxfile05 = "mvtx5.list",
0100     const std::string &input_tpcfile00_0 = "tpc00_0.list",
0101     const std::string &input_tpcfile00_1 = "tpc00_1.list",
0102     const std::string &input_tpcfile01_0 = "tpc01_0.list",
0103     const std::string &input_tpcfile01_1 = "tpc01_1.list",
0104     const std::string &input_tpcfile02_0 = "tpc02_0.list",
0105     const std::string &input_tpcfile02_1 = "tpc02_1.list",
0106     const std::string &input_tpcfile03_0 = "tpc03_0.list",
0107     const std::string &input_tpcfile03_1 = "tpc03_1.list",
0108     const std::string &input_tpcfile04_0 = "tpc04_0.list",
0109     const std::string &input_tpcfile04_1 = "tpc04_1.list",
0110     const std::string &input_tpcfile05_0 = "tpc05_0.list",
0111     const std::string &input_tpcfile05_1 = "tpc05_1.list",
0112     const std::string &input_tpcfile06_0 = "tpc06_0.list",
0113     const std::string &input_tpcfile06_1 = "tpc06_1.list",
0114     const std::string &input_tpcfile07_0 = "tpc07_0.list",
0115     const std::string &input_tpcfile07_1 = "tpc07_1.list",
0116     const std::string &input_tpcfile08_0 = "tpc08_0.list",
0117     const std::string &input_tpcfile08_1 = "tpc08_1.list",
0118     const std::string &input_tpcfile09_0 = "tpc09_0.list",
0119     const std::string &input_tpcfile09_1 = "tpc09_1.list",
0120     const std::string &input_tpcfile10_0 = "tpc10_0.list",
0121     const std::string &input_tpcfile10_1 = "tpc10_1.list",
0122     const std::string &input_tpcfile11_0 = "tpc11_0.list",
0123     const std::string &input_tpcfile11_1 = "tpc11_1.list",
0124     const std::string &input_tpcfile12_0 = "tpc12_0.list",
0125     const std::string &input_tpcfile12_1 = "tpc12_1.list",
0126     const std::string &input_tpcfile13_0 = "tpc13_0.list",
0127     const std::string &input_tpcfile13_1 = "tpc13_1.list",
0128     const std::string &input_tpcfile14_0 = "tpc14_0.list",
0129     const std::string &input_tpcfile14_1 = "tpc14_1.list",
0130     const std::string &input_tpcfile15_0 = "tpc15_0.list",
0131     const std::string &input_tpcfile15_1 = "tpc15_1.list",
0132     const std::string &input_tpcfile16_0 = "tpc16_0.list",
0133     const std::string &input_tpcfile16_1 = "tpc16_1.list",
0134     const std::string &input_tpcfile17_0 = "tpc17_0.list",
0135     const std::string &input_tpcfile17_1 = "tpc17_1.list",
0136     const std::string &input_tpcfile18_0 = "tpc18_0.list",
0137     const std::string &input_tpcfile18_1 = "tpc18_1.list",
0138     const std::string &input_tpcfile19_0 = "tpc19_0.list",
0139     const std::string &input_tpcfile19_1 = "tpc19_1.list",
0140     const std::string &input_tpcfile20_0 = "tpc20_0.list",
0141     const std::string &input_tpcfile20_1 = "tpc20_1.list",
0142     const std::string &input_tpcfile21_0 = "tpc21_0.list",
0143     const std::string &input_tpcfile21_1 = "tpc21_1.list",
0144     const std::string &input_tpcfile22_0 = "tpc22_0.list",
0145     const std::string &input_tpcfile22_1 = "tpc22_1.list",
0146     const std::string &input_tpcfile23_0 = "tpc23_0.list",
0147     const std::string &input_tpcfile23_1 = "tpc23_1.list",
0148     const std::string &input_tpotfile = "tpot.list")
0149 {
0150   std::vector<std::string> gl1_infile;
0151   gl1_infile.push_back(input_gl1file);
0152 
0153   // MVTX
0154   std::vector<std::string> mvtx_infile;
0155   mvtx_infile.push_back(input_mvtxfile00);
0156   mvtx_infile.push_back(input_mvtxfile01);
0157   mvtx_infile.push_back(input_mvtxfile02);
0158   mvtx_infile.push_back(input_mvtxfile03);
0159   mvtx_infile.push_back(input_mvtxfile04);
0160   mvtx_infile.push_back(input_mvtxfile05);
0161 
0162   // INTT
0163   std::vector<std::string> intt_infile;
0164   intt_infile.push_back(input_inttfile00);
0165   intt_infile.push_back(input_inttfile01);
0166   intt_infile.push_back(input_inttfile02);
0167   intt_infile.push_back(input_inttfile03);
0168   intt_infile.push_back(input_inttfile04);
0169   intt_infile.push_back(input_inttfile05);
0170   intt_infile.push_back(input_inttfile06);
0171   intt_infile.push_back(input_inttfile07);
0172 
0173   std::vector<std::string> tpc_infile;
0174   tpc_infile.push_back(input_tpcfile00_0);
0175   tpc_infile.push_back(input_tpcfile01_0);
0176   tpc_infile.push_back(input_tpcfile02_0);
0177   tpc_infile.push_back(input_tpcfile03_0);
0178   tpc_infile.push_back(input_tpcfile04_0);
0179   tpc_infile.push_back(input_tpcfile05_0);
0180   tpc_infile.push_back(input_tpcfile06_0);
0181   tpc_infile.push_back(input_tpcfile07_0);
0182   tpc_infile.push_back(input_tpcfile08_0);
0183   tpc_infile.push_back(input_tpcfile09_0);
0184   tpc_infile.push_back(input_tpcfile10_0);
0185   tpc_infile.push_back(input_tpcfile11_0);
0186   tpc_infile.push_back(input_tpcfile12_0);
0187   tpc_infile.push_back(input_tpcfile13_0);
0188   tpc_infile.push_back(input_tpcfile14_0);
0189   tpc_infile.push_back(input_tpcfile15_0);
0190   tpc_infile.push_back(input_tpcfile16_0);
0191   tpc_infile.push_back(input_tpcfile17_0);
0192   tpc_infile.push_back(input_tpcfile18_0);
0193   tpc_infile.push_back(input_tpcfile19_0);
0194   tpc_infile.push_back(input_tpcfile20_0);
0195   tpc_infile.push_back(input_tpcfile21_0);
0196   tpc_infile.push_back(input_tpcfile22_0);
0197   tpc_infile.push_back(input_tpcfile23_0);
0198   tpc_infile.push_back(input_tpcfile00_1);
0199   tpc_infile.push_back(input_tpcfile01_1);
0200   tpc_infile.push_back(input_tpcfile02_1);
0201   tpc_infile.push_back(input_tpcfile03_1);
0202   tpc_infile.push_back(input_tpcfile04_1);
0203   tpc_infile.push_back(input_tpcfile05_1);
0204   tpc_infile.push_back(input_tpcfile06_1);
0205   tpc_infile.push_back(input_tpcfile07_1);
0206   tpc_infile.push_back(input_tpcfile08_1);
0207   tpc_infile.push_back(input_tpcfile09_1);
0208   tpc_infile.push_back(input_tpcfile10_1);
0209   tpc_infile.push_back(input_tpcfile11_1);
0210   tpc_infile.push_back(input_tpcfile12_1);
0211   tpc_infile.push_back(input_tpcfile13_1);
0212   tpc_infile.push_back(input_tpcfile14_1);
0213   tpc_infile.push_back(input_tpcfile15_1);
0214   tpc_infile.push_back(input_tpcfile16_1);
0215   tpc_infile.push_back(input_tpcfile17_1);
0216   tpc_infile.push_back(input_tpcfile18_1);
0217   tpc_infile.push_back(input_tpcfile19_1);
0218   tpc_infile.push_back(input_tpcfile20_1);
0219   tpc_infile.push_back(input_tpcfile21_1);
0220   tpc_infile.push_back(input_tpcfile22_1);
0221   tpc_infile.push_back(input_tpcfile23_1);
0222 
0223   // TPOT
0224   std::vector<std::string> tpot_infile;
0225   tpot_infile.push_back(input_tpotfile);
0226 
0227   int runnumber = -99999;
0228   if (!gl1_infile.empty())
0229   {
0230     runnumber = getrunnumber(gl1_infile[0]);
0231   }
0232   else if (!mvtx_infile.empty())
0233   {
0234     runnumber = getrunnumber(mvtx_infile[0]);
0235   }
0236   else if (!intt_infile.empty())
0237   {
0238     runnumber = getrunnumber(intt_infile[0]);
0239   }
0240   else if (!tpc_infile.empty())
0241   {
0242     runnumber = getrunnumber(tpc_infile[0]);
0243   }
0244   else if (!tpot_infile.empty())
0245   {
0246     runnumber = getrunnumber(tpot_infile[0]);
0247   }
0248   if (runnumber == -99999)
0249   {
0250     std::cout << "could not extract run number from input files (all lists empty?)"
0251      << std::endl;
0252     gSystem->Exit(1);
0253   }
0254   auto *se = Fun4AllServer::instance();
0255   se->Verbosity(2);
0256   auto *rc = recoConsts::instance();
0257   rc->set_IntFlag("RUNNUMBER", runnumber);
0258 
0259   Enable::CDB = true;
0260   rc->set_StringFlag("CDB_GLOBALTAG", "newcdbtag");
0261 
0262   rc->set_uint64Flag("TIMESTAMP", runnumber);
0263 
0264   //! flags to set
0265   Enable::QA = true;
0266   TRACKING::tpc_zero_supp = true;
0267   G4TRACKING::convert_seeds_to_svtxtracks = false;
0268 
0269   Enable::MVTX_APPLYMISALIGNMENT = true;
0270   ACTSGEOM::mvtx_applymisalignment = Enable::MVTX_APPLYMISALIGNMENT;
0271   TpcReadoutInit(runnumber);
0272   std::cout << " run: " << runnumber
0273             << " samples: " << TRACKING::reco_tpc_maxtime_sample
0274             << " pre: " << TRACKING::reco_tpc_time_presample
0275             << " vdrift: " << G4TPC::tpc_drift_velocity_reco
0276             << std::endl;
0277 
0278   CDBInterface::instance()->Verbosity(1);
0279   std::string geofile = CDBInterface::instance()->getUrl("Tracking_Geometry");
0280   std::cout << "CDB tracking geometry file " << geofile << std::endl;
0281   Fun4AllRunNodeInputManager *ingeo = new Fun4AllRunNodeInputManager("GeoIn");
0282   ingeo->AddFile(geofile);
0283   se->registerInputManager(ingeo);
0284 
0285   // can use for zero field
0286   // double fieldstrength = 0.01;
0287   // G4MAGNET::magfield_tracking = "0.01";
0288   double fieldstrength = std::numeric_limits<double>::quiet_NaN();
0289   bool ConstField = isConstantField(G4MAGNET::magfield_tracking, fieldstrength);
0290 
0291   if (ConstField && fieldstrength < 0.1)
0292   {
0293     G4MAGNET::magfield = "0.01";
0294     G4MAGNET::magfield_rescale = 1;
0295   }
0296 
0297   TrackingInit();
0298 
0299   int i = 0;
0300   int NumInputs = 0;
0301   Fun4AllStreamingInputManager *in = new Fun4AllStreamingInputManager("Comb");
0302 
0303   for (const auto& iter : gl1_infile)
0304   {
0305     if (isGood(iter))
0306     {
0307       SingleGl1PoolInput *gl1_sngl = new SingleGl1PoolInput("GL1_" + std::to_string(i));
0308       //    gl1_sngl->Verbosity(3);
0309       gl1_sngl->AddListFile(iter);
0310       in->registerStreamingInput(gl1_sngl, InputManagerType::GL1);
0311       i++;
0312     }
0313   }
0314   NumInputs += i;
0315 
0316   i = 0;
0317   for (const auto& iter : intt_infile)
0318   {
0319     if (isGood(iter))
0320     {
0321       SingleInttPoolInput *intt_sngl = new SingleInttPoolInput("INTT_" + std::to_string(i));
0322       //    intt_sngl->Verbosity(3);
0323       intt_sngl->AddListFile(iter);
0324       in->registerStreamingInput(intt_sngl, InputManagerType::INTT);
0325       i++;
0326     }
0327   }
0328   NumInputs += i;
0329 
0330   i = 0;
0331   for (const auto& iter : mvtx_infile)
0332   {
0333     if (isGood(iter))
0334     {
0335       SingleMvtxPoolInput *mvtx_sngl = new SingleMvtxPoolInput("MVTX_" + std::to_string(i));
0336       //    mvtx_sngl->Verbosity(3);
0337       mvtx_sngl->AddListFile(iter);
0338       in->registerStreamingInput(mvtx_sngl, InputManagerType::MVTX);
0339       i++;
0340     }
0341   }
0342   NumInputs += i;
0343 
0344   i = 0;
0345 
0346   for (const auto& iter : tpc_infile)
0347   {
0348     if (isGood(iter))
0349     {
0350       SingleTpcTimeFrameInput *tpc_sngl = new SingleTpcTimeFrameInput("TPC_" + std::to_string(i));
0351       //    tpc_sngl->Verbosity(2);
0352       //   tpc_sngl->DryRun();
0353       tpc_sngl->setHitContainerName("TPCRAWHIT");
0354       tpc_sngl->AddListFile(iter);
0355       in->registerStreamingInput(tpc_sngl, InputManagerType::TPC);
0356       i++;
0357     }
0358   }
0359   NumInputs += i;
0360   i = 0;
0361   for (const auto& iter : tpot_infile)
0362   {
0363     if (isGood(iter))
0364     {
0365       SingleMicromegasPoolInput *mm_sngl = new SingleMicromegasPoolInput("MICROMEGAS_" + std::to_string(i));
0366       //   sngl->Verbosity(3);
0367       mm_sngl->SetBcoRange(10);
0368       mm_sngl->SetNegativeBco(2);
0369       mm_sngl->SetBcoPoolSize(50);
0370       mm_sngl->AddListFile(iter);
0371       in->registerStreamingInput(mm_sngl, InputManagerType::MICROMEGAS);
0372       i++;
0373     }
0374   }
0375   NumInputs += i;
0376 
0377   // if there is no input manager this macro will still run - so just quit here
0378   if (NumInputs == 0)
0379   {
0380     std::cout << "no file lists no input manager registered, quitting" << std::endl;
0381     gSystem->Exit(1);
0382   }
0383   se->registerInputManager(in);
0384 
0385   SyncReco *sync = new SyncReco();
0386   se->registerSubsystem(sync);
0387 
0388   HeadReco *head = new HeadReco();
0389   se->registerSubsystem(head);
0390 
0391   FlagHandler *flag = new FlagHandler();
0392   se->registerSubsystem(flag);
0393 
0394   Mvtx_HitUnpacking();
0395   Intt_HitUnpacking();
0396   Tpc_HitUnpacking();
0397   Micromegas_HitUnpacking();
0398 
0399   MvtxClusterizer *mvtxclusterizer = new MvtxClusterizer("MvtxClusterizer");
0400 
0401   se->registerSubsystem(mvtxclusterizer);
0402 
0403   Intt_Clustering();
0404 
0405   Tpc_LaserEventIdentifying();
0406 
0407   auto *tpcclusterizer = new TpcClusterizer;
0408   tpcclusterizer->Verbosity(0);
0409   tpcclusterizer->set_do_hit_association(G4TPC::DO_HIT_ASSOCIATION);
0410   tpcclusterizer->set_rawdata_reco();
0411   tpcclusterizer->set_reject_event(G4TPC::REJECT_LASER_EVENTS);
0412   se->registerSubsystem(tpcclusterizer);
0413 
0414   Micromegas_Clustering();
0415 
0416   auto *silicon_Seeding = new PHActsSiliconSeeding;
0417   silicon_Seeding->Verbosity(0);
0418   // these get us to about 83% INTT > 1
0419   silicon_Seeding->setinttRPhiSearchWindow(0.4);
0420   silicon_Seeding->setinttZSearchWindow(2.0);
0421   silicon_Seeding->setStrobeRange(-5, 5);
0422   silicon_Seeding->seedAnalysis(false);
0423   se->registerSubsystem(silicon_Seeding);
0424 
0425   auto *merger = new PHSiliconSeedMerger;
0426   merger->Verbosity(0);
0427   se->registerSubsystem(merger);
0428 
0429   /*
0430    * Tpc Seeding
0431    */
0432   auto *seeder = new PHCASeeding("PHCASeeding");
0433   if (ConstField)
0434   {
0435     seeder->useConstBField(true);
0436     seeder->constBField(fieldstrength);
0437   }
0438   else
0439   {
0440     seeder->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0441     seeder->useConstBField(false);
0442     seeder->magFieldFile(G4MAGNET::magfield_tracking);  // to get charge sign right
0443   }
0444   seeder->Verbosity(0);
0445   seeder->SetLayerRange(7, 55);
0446   seeder->SetSearchWindow(2., 0.05);           // z-width and phi-width, default in macro at 1.5 and 0.05
0447   seeder->SetClusAdd_delta_window(3.0, 0.06);  //  (0.5, 0.005) are default; sdzdr_cutoff, d2/dr2(phi)_cutoff
0448   // seeder->SetNClustersPerSeedRange(4,60); // default is 6, 6
0449   seeder->SetMinHitsPerCluster(0);
0450   seeder->SetMinClustersPerTrack(3);
0451   seeder->useFixedClusterError(true);
0452   seeder->set_pp_mode(true);
0453   se->registerSubsystem(seeder);
0454 
0455   // expand stubs in the TPC using simple kalman filter
0456   auto *cprop = new PHSimpleKFProp("PHSimpleKFProp");
0457   cprop->set_field_dir(G4MAGNET::magfield_rescale);
0458   if (ConstField)
0459   {
0460     cprop->useConstBField(true);
0461     cprop->setConstBField(fieldstrength);
0462   }
0463   else
0464   {
0465     cprop->magFieldFile(G4MAGNET::magfield_tracking);
0466     cprop->set_field_dir(-1 * G4MAGNET::magfield_rescale);
0467   }
0468   cprop->useFixedClusterError(true);
0469   cprop->set_max_window(5.);
0470   cprop->Verbosity(0);
0471   cprop->set_pp_mode(true);
0472   se->registerSubsystem(cprop);
0473 
0474   // Always apply preliminary distortion corrections to TPC clusters before silicon matching
0475   // and refit the trackseeds. Replace KFProp fits with the new fit parameters in the TPC seeds.
0476   auto *prelim_distcorr = new PrelimDistortionCorrection;
0477   prelim_distcorr->set_pp_mode(true);
0478   prelim_distcorr->Verbosity(0);
0479   se->registerSubsystem(prelim_distcorr);
0480 
0481   /*
0482    * Track Matching between silicon and TPC
0483    */
0484   // The normal silicon association methods
0485   // Match the TPC track stubs from the CA seeder to silicon track stubs from PHSiliconTruthTrackSeeding
0486   auto *silicon_match = new PHSiliconTpcTrackMatching;
0487   silicon_match->Verbosity(0);
0488   silicon_match->set_pp_mode(TRACKING::pp_mode);
0489   if (G4TPC::ENABLE_AVERAGE_CORRECTIONS)
0490   {
0491     // for general tracking
0492     // Eta/Phi window is determined by 3 sigma window
0493     // X/Y/Z window is determined by 4 sigma window
0494     silicon_match->window_deta.set_posQoverpT_maxabs({-0.014, 0.0331, 0.48});
0495     silicon_match->window_deta.set_negQoverpT_maxabs({-0.006, 0.0235, 0.52});
0496     silicon_match->set_deltaeta_min(0.03);
0497     silicon_match->window_dphi.set_QoverpT_range({-0.15, 0, 0}, {0.15, 0, 0});
0498     silicon_match->window_dx.set_QoverpT_maxabs({3.0, 0, 0});
0499     silicon_match->window_dy.set_QoverpT_maxabs({3.0, 0, 0});
0500     silicon_match->window_dz.set_posQoverpT_maxabs({1.138, 0.3919, 0.84});
0501     silicon_match->window_dz.set_negQoverpT_maxabs({0.719, 0.6485, 0.65});
0502     silicon_match->set_crossing_deltaz_max(30);
0503     silicon_match->set_crossing_deltaz_min(2);
0504 
0505     // for distortion correction using SI-TPOT fit and track pT>0.5
0506     if (G4TRACKING::SC_CALIBMODE)
0507     {
0508       silicon_match->window_deta.set_posQoverpT_maxabs({0.016, 0.0060, 1.13});
0509       silicon_match->window_deta.set_negQoverpT_maxabs({0.022, 0.0022, 1.44});
0510       silicon_match->set_deltaeta_min(0.03);
0511       silicon_match->window_dphi.set_QoverpT_range({-0.15, 0, 0}, {0.09, 0, 0});
0512       silicon_match->window_dx.set_QoverpT_maxabs({2.0, 0, 0});
0513       silicon_match->window_dy.set_QoverpT_maxabs({1.5, 0, 0});
0514       silicon_match->window_dz.set_posQoverpT_maxabs({1.213, 0.0211, 2.09});
0515       silicon_match->window_dz.set_negQoverpT_maxabs({1.307, 0.0001, 4.52});
0516       silicon_match->set_crossing_deltaz_min(1.2);
0517     }
0518   }
0519   se->registerSubsystem(silicon_match);
0520 
0521   // Match TPC track stubs from CA seeder to clusters in the micromegas layers
0522   auto *mm_match = new PHMicromegasTpcTrackMatching;
0523   mm_match->Verbosity(0);
0524   mm_match->set_rphi_search_window_lyr1(3.);
0525   mm_match->set_rphi_search_window_lyr2(15.0);
0526   mm_match->set_z_search_window_lyr1(30.0);
0527   mm_match->set_z_search_window_lyr2(3.);
0528 
0529   mm_match->set_min_tpc_layer(38);             // layer in TPC to start projection fit
0530   mm_match->set_test_windows_printout(false);  // used for tuning search windows only
0531   se->registerSubsystem(mm_match);
0532 
0533   /*
0534    * End Track Seeding
0535    */
0536 
0537   /*
0538    * Either converts seeds to tracks with a straight line/helix fit
0539    * or run the full Acts track kalman filter fit
0540    */
0541   if (G4TRACKING::convert_seeds_to_svtxtracks)
0542   {
0543     auto *converter = new TrackSeedTrackMapConverter;
0544     // Default set to full SvtxTrackSeeds. Can be set to
0545     // SiliconTrackSeedContainer or TpcTrackSeedContainer
0546     converter->setTrackSeedName("SvtxTrackSeedContainer");
0547     converter->setFieldMap(G4MAGNET::magfield_tracking);
0548     converter->Verbosity(0);
0549     se->registerSubsystem(converter);
0550   }
0551   else
0552   {
0553     auto *deltazcorr = new PHTpcDeltaZCorrection;
0554     deltazcorr->Verbosity(0);
0555     se->registerSubsystem(deltazcorr);
0556 
0557     // perform final track fit with ACTS
0558     auto *actsFit = new PHActsTrkFitter;
0559     actsFit->Verbosity(0);
0560     actsFit->commissioning(G4TRACKING::use_alignment);
0561     // in calibration mode, fit only Silicons and Micromegas hits
0562     actsFit->fitSiliconMMs(G4TRACKING::SC_CALIBMODE);
0563     actsFit->setUseMicromegas(G4TRACKING::SC_USE_MICROMEGAS);
0564     actsFit->set_pp_mode(TRACKING::pp_mode);
0565     actsFit->set_use_clustermover(true);  // default is true for now
0566     actsFit->useActsEvaluator(false);
0567     actsFit->useOutlierFinder(false);
0568     actsFit->setFieldMap(G4MAGNET::magfield_tracking);
0569     se->registerSubsystem(actsFit);
0570 
0571     auto *cleaner = new PHTrackCleaner();
0572     cleaner->Verbosity(0);
0573     cleaner->set_pp_mode(TRACKING::pp_mode);
0574     // se->registerSubsystem(cleaner);
0575   }
0576 
0577   auto *finder = new PHSimpleVertexFinder;
0578   finder->Verbosity(0);
0579   finder->setDcaCut(0.5);
0580   finder->setTrackPtCut(-99999.);
0581   finder->setBeamLineCut(1);
0582   finder->setTrackQualityCut(1000000000);
0583   finder->setNmvtxRequired(3);
0584   finder->setOutlierPairCut(0.1);
0585   se->registerSubsystem(finder);
0586 
0587   std::string residstring = outfilename + "_resid.root";
0588 
0589   auto *resid = new TrackResiduals("TrackResiduals");
0590   resid->Verbosity(0);
0591   resid->outfileName(residstring);
0592   resid->alignment(false);
0593   resid->clusterTree();
0594   // resid->failedTree();
0595   // resid->hitTree();
0596   // resid->noEventTree(); // not implemented anymore?
0597   resid->convertSeeds(G4TRACKING::convert_seeds_to_svtxtracks);
0598 
0599   if (ConstField && fieldstrength < 0.1)
0600   {
0601     resid->zeroField();
0602   }
0603   se->registerSubsystem(resid);
0604 
0605   // Fun4AllOutputManager *out = new Fun4AllDstOutputManager("out", "/sphenix/tg/tg01/hf/jdosbo/tracking_development/onlineoffline/hitsets.root");
0606   // se->registerOutputManager(out);
0607 
0608   if (Enable::QA)
0609   {
0610     se->registerSubsystem(new MvtxRawHitQA);
0611     se->registerSubsystem(new InttRawHitQA);
0612     se->registerSubsystem(new TpcRawHitQA);
0613     se->registerSubsystem(new MvtxClusterQA);
0614     se->registerSubsystem(new InttClusterQA);
0615     se->registerSubsystem(new TpcClusterQA);
0616     se->registerSubsystem(new MicromegasClusterQA);
0617     se->registerSubsystem(new SiliconSeedsQA);
0618     se->registerSubsystem(new TpcSeedsQA);
0619     se->registerSubsystem(new TpcSiliconQA);
0620   }
0621 
0622   se->run(nEvents);
0623   se->End();
0624   se->PrintTimer();
0625 
0626   if (Enable::QA)
0627   {
0628     TString qaname = outfilename + runnumber + "_qa.root";
0629     std::string qaOutputFileName = outfilename + std::to_string(runnumber) + "_qa.root";
0630     QAHistManagerDef::saveQARootFile(qaOutputFileName);
0631   }
0632   CDBInterface::instance()->Print();
0633   delete se;
0634   std::cout << "Finished" << std::endl;
0635   gSystem->Exit(0);
0636 }
0637 
0638 bool isGood(const std::string &infile)
0639 {
0640   std::ifstream intest;
0641   intest.open(infile);
0642   bool goodfile = false;
0643   if (intest.is_open())
0644   {
0645     if (intest.peek() != std::ifstream::traits_type::eof())  // is it non zero?
0646     {
0647       goodfile = true;
0648     }
0649     intest.close();
0650   }
0651   return goodfile;
0652 }
0653 
0654 int getrunnumber(const std::string &listfile)
0655 {
0656   if (! isGood(listfile))
0657   {
0658     std::cout << "listfile " << listfile << " is bad" << std::endl;
0659     gSystem->Exit(1);
0660   }
0661   std::ifstream ifs(listfile);
0662   std::string filepath;
0663   std::getline(ifs, filepath);
0664 
0665   std::pair<int, int> runseg = Fun4AllUtils::GetRunSegment(filepath);
0666   int runnumber = runseg.first;
0667 //  int segment = abs(runseg.second);
0668   return runnumber;
0669 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!