macros/data/Fun4All_HF_TPConly.C

0001 #include <fun4all/Fun4AllUtils.h>
0002 #include <G4_ActsGeom.C>
0003 #include <G4_Global.C>
0004 #include <G4_Magnet.C>
0005 #include <G4_Mbd.C>
0006 #include <GlobalVariables.C>
0007 #include <Trkr_Clustering.C>
0008 #include <Trkr_LaserClustering.C>
0009 #include <Trkr_Reco.C>
0010 #include <Trkr_RecoInit.C>
0011 #include <Trkr_TpcReadoutInit.C>
0012 #include <QA.C>
0013 #include <Calo_Calib.C>
0014
0015 #include <ffamodules/CDBInterface.h>
0016
0017 #include <fun4all/Fun4AllDstInputManager.h>
0018 #include <fun4all/Fun4AllDstOutputManager.h>
0019 #include <fun4all/Fun4AllInputManager.h>
0020 #include <fun4all/Fun4AllOutputManager.h>
0021 #include <fun4all/Fun4AllRunNodeInputManager.h>
0022 #include <fun4all/Fun4AllServer.h>
0023
0024 #include <phool/recoConsts.h>
0025
0026 #include <mvtxrawhitqa/MvtxRawHitQA.h>
0027 #include <inttrawhitqa/InttRawHitQA.h>
0028 #include <tpcqa/TpcRawHitQA.h>
0029 #include <trackingqa/InttClusterQA.h>
0030 #include <trackingqa/MicromegasClusterQA.h>
0031 #include <trackingqa/MvtxClusterQA.h>
0032 #include <trackingqa/TpcClusterQA.h>
0033 #include <trackingqa/SiliconSeedsQA.h>
0034 #include <trackingqa/TpcSeedsQA.h>
0035 #include <trackingqa/TrackFittingQA.h>
0036 #include <trackingqa/TpcSiliconQA.h>
0037 #include <trackingqa/VertexQA.h>
0038 #include <trackreco/PHActsTrackProjection.h>
0039
0040 #include <caloreco/CaloGeomMapping.h>
0041 #include <caloreco/RawClusterBuilderTemplate.h>
0042
0043 #include <kfparticle_sphenix/KFParticle_sPHENIX.h>
0044
0045 #include <stdio.h>
0046 //#include <nnscalemap/nnscalemap.h>
0047 #include "HF_selections_TPConly.C"
0048
0049
0050 using namespace HeavyFlavorReco;
0051
0052 R__LOAD_LIBRARY(libfun4all.so)
0053 //R__LOAD_LIBRARY(libNNScaleMap.so)
0054 R__LOAD_LIBRARY(libffamodules.so)
0055 R__LOAD_LIBRARY(libphool.so)
0056 R__LOAD_LIBRARY(libcdbobjects.so)
0057
0058 bool IsListFile = false;
0059 bool IsDstFile = false;
0060 bool DoUnpacking = false;
0061 bool DoClustering = false;
0062 bool DoSeeding = false;
0063 bool DoFitting = false;
0064
0065 void CheckDstType(const std::string inputDST)
0066 {
0067     if (inputDST.find("DST_STREAMING_EVENT") != std::string::npos)
0068     {
0069         DoUnpacking = true;
0070         DoClustering = true;
0071         DoSeeding = true;
0072         DoFitting = true;
0073     }
0074     else if (inputDST.find("DST_TRKR_CLUSTER") != std::string::npos)
0075     {
0076         DoUnpacking = false;
0077         DoClustering = false;
0078         DoSeeding = true;
0079         DoFitting = true;
0080     }
0081     else if (inputDST.find("DST_TRKR_SEED") != std::string::npos)
0082     {
0083         DoUnpacking = false;
0084         DoClustering = false;
0085         DoSeeding = false;
0086         DoFitting = true;
0087     }
0088     else if (inputDST.find("DST_TRKR_TRACKS") != std::string::npos)
0089     {
0090         DoUnpacking = false;
0091         DoClustering = false;
0092         DoSeeding = false;
0093         DoFitting = false;
0094     }
0095     return;
0096 }
0097
0098 std::string GetRunInfo(const std::string &filename)
0099 {
0100     const std::vector<std::string> runspecies_patterns = {"run2pp", "run3pp", "run3auau", "run3oo"};
0101
0102     boost::char_separator<char> sep("_");
0103     boost::tokenizer<boost::char_separator<char>> tok(filename, sep);
0104
0105     for (const auto &t : tok)
0106     {
0107         for (const auto &pattern : runspecies_patterns)
0108         {
0109             if (t == pattern)
0110             {
0111                 return t;
0112             }
0113         }
0114     }
0115
0116     std::cout << "Cannot extract run info from filename!" << std::endl;
0117     return "";
0118 }
0119
0120 void Fun4All_HF_TPConly(const int nEvents = -1,
0121                 const std::string inputDST = "/sphenix/lustre01/sphnxpro/production/run3pp/physics/ana538_2025p011_v001/DST_TRKR_SEED/run_00079500_00079600/DST_TRKR_SEED_run3pp_ana538_2025p011_v001-00079516-01995.root", //
0122                 const std::string outDir = "./",
0123                 const int nSkip = 0,                                                                                                                                                                                                //
0124                 const bool convertSeeds = true)                                                                                                                                                                                    //
0125 {
0126     auto se = Fun4AllServer::instance();
0127     se->Verbosity(1);
0128     // se->VerbosityDownscale(100);
0129     auto rc = recoConsts::instance();
0130
0131     if (inputDST.find(".list") != std::string::npos)
0132     {
0133         IsListFile = true;
0134         std::cout << "Input is list file" << std::endl;
0135     }
0136     if (inputDST.find(".root") != std::string::npos)
0137     {
0138         IsDstFile = true;
0139         std::cout << "Input is dst file" << std::endl;
0140     }
0141     if (!IsListFile && !IsDstFile)
0142     {
0143         std::cout << "Check your input! Exit" << std::endl;
0144         gSystem->Exit(0);
0145     }
0146
0147     int runnumber = std::numeric_limits<int>::quiet_NaN();
0148     int segment = std::numeric_limits<int>::quiet_NaN();
0149     std::string runspecies = "";
0150
0151     if (IsListFile)
0152     {
0153         std::ifstream ifs(inputDST);
0154         std::string filepath;
0155         int i = 0;
0156         while (std::getline(ifs, filepath))
0157         {
0158             std::cout << "Adding DST with filepath: " << filepath << std::endl;
0159             if (i == 0)
0160             {
0161                 std::pair<int, int> runseg = Fun4AllUtils::GetRunSegment(filepath);
0162                 runnumber = runseg.first;
0163                 segment = runseg.second;
0164                 runspecies = GetRunInfo(filepath);
0165                 CheckDstType(filepath);
0166             }
0167             std::string inputname = "InputManager" + std::to_string(i);
0168             auto *hitsin = new Fun4AllDstInputManager(inputname);
0169             hitsin->fileopen(filepath);
0170             se->registerInputManager(hitsin);
0171             i++;
0172         }
0173     }
0174
0175     if (IsDstFile)
0176     {
0177         std::pair<int, int> runseg = Fun4AllUtils::GetRunSegment(inputDST);
0178         runnumber = runseg.first;
0179         segment = runseg.second;
0180         runspecies = GetRunInfo(inputDST);
0181
0182         CheckDstType(inputDST);
0183
0184         auto *hitsin = new Fun4AllDstInputManager("InputManager");
0185         hitsin->fileopen(inputDST);
0186         se->registerInputManager(hitsin);
0187     }
0188
0189     rc->set_IntFlag("RUNNUMBER", runnumber);
0190     rc->set_IntFlag("RUNSEGMENT", segment);
0191     rc->set_uint64Flag("TIMESTAMP", runnumber);
0192
0193     Enable::CDB = true;
0194     rc->set_StringFlag("CDB_GLOBALTAG", "newcdbtag");
0195     rc->set_uint64Flag("TIMESTAMP", runnumber);
0196
0197     std::stringstream nice_runnumber;
0198     nice_runnumber << std::setw(8) << std::setfill('0') << to_string(runnumber);
0199
0200     int rounded_up = 100 * (std::ceil((float)(runnumber+0.1) / 100));
0201     std::stringstream nice_rounded_up;
0202     nice_rounded_up << std::setw(8) << std::setfill('0') << to_string(rounded_up);
0203
0204     int rounded_down = 100 * (std::floor((float)(runnumber+0.1) / 100));
0205     std::stringstream nice_rounded_down;
0206     nice_rounded_down << std::setw(8) << std::setfill('0') << to_string(rounded_down);
0207
0208     std::stringstream nice_segment;
0209     nice_segment << std::setw(5) << std::setfill('0') << to_string(segment);
0210
0211     std::stringstream nice_skip;
0212     nice_skip << std::setw(5) << std::setfill('0') << to_string(nSkip);
0213
0214
0215     if (!DoSeeding && (get_dEdx_info || get_detector_info))
0216     {
0217         //std::string clus_anacdbver = "ana534_2025p009_v001";
0218         //std::string clus_anacdbver = "new_newcdbtag_v001";
0219         //std::string clus_anacdbver = "ana532_2025p009_v001";
0220         std::string clus_anacdbver = "ana538_2025p011_v001";
0221
0222         std::string clus_file = "";
0223         if (runspecies == "run3pp" || runspecies == "run3auau" || runspecies == "run3oo")
0224         {
0225             clus_file = "/sphenix/lustre01/sphnxpro/production/" + runspecies + "/physics/" + clus_anacdbver + "/DST_TRKR_CLUSTER/run_" + nice_rounded_down.str() + "_" + nice_rounded_up.str() + "/DST_TRKR_CLUSTER_" + runspecies + "_" + clus_anacdbver + "-" + nice_runnumber.str() + "-" + nice_segment.str() + ".root";
0226         }
0227         if (runspecies == "run2pp")
0228         {
0229             clus_file = "/sphenix/lustre01/sphnxpro/production/" + runspecies + "/physics/" + clus_anacdbver + "/DST_TRKR_CLUSTER/run_" + nice_rounded_down.str() + "_" + nice_rounded_up.str() + "/dst/DST_TRKR_CLUSTER_" + runspecies + "_" + clus_anacdbver + "-" + nice_runnumber.str() + "-" + nice_segment.str() + ".root";
0230         }
0231         std::cout << "Input cluster DST: " << clus_file << std::endl;
0232
0233         auto hitsinclus = new Fun4AllDstInputManager("ClusterInputManager");
0234         hitsinclus->fileopen(clus_file);
0235         se->registerInputManager(hitsinclus);
0236     }
0237
0238
0239     /*
0240     if (get_calo_info)
0241     {
0242         int ratio_nCalo_over_nTrkr = 10; // number of events per calo DST divided by number of events per tracking DST
0243         std::stringstream nice_segment_calo;
0244         nice_segment_calo << std::setw(5) << std::setfill('0') << to_string(segment / ratio_nCalo_over_nTrkr);
0245
0246         std::string calo_anacdbver = "ana509_2024p022_v001";
0247         std::string calo_file = "";
0248         calo_file = "/sphenix/lustre01/sphnxpro/production2/" + runspecies + "/physics/calofitting/" + calo_anacdbver + "/run_" + nice_rounded_down.str() + "_" + nice_rounded_up.str() + "/DST_CALOFITTING_" + runspecies + "_" + calo_anacdbver + "-" + nice_runnumber.str() + "-" + nice_segment_calo.str() + ".root";
0249         std::cout << "Input calo DST: " << calo_file << std::endl;
0250
0251         auto hitsincalo = new Fun4AllDstInputManager("CaloInputManager");
0252         hitsincalo->fileopen(calo_file);
0253         se->registerInputManager(hitsincalo);
0254     }
0255     */
0256
0257     std::string geofile = CDBInterface::instance()->getUrl("Tracking_Geometry");
0258
0259     Fun4AllRunNodeInputManager *ingeo = new Fun4AllRunNodeInputManager("GeoIn");
0260     ingeo->AddFile(geofile);
0261     se->registerInputManager(ingeo);
0262
0263     // set flags
0264     //TRACKING::pp_mode = true;
0265     TRACKING::streaming_mode = true;
0266
0267     Enable::MVTX_APPLYMISALIGNMENT = true;
0268     ACTSGEOM::mvtx_applymisalignment = Enable::MVTX_APPLYMISALIGNMENT;
0269
0270     // distortion calibration mode
0271     /*
0272      * set to true to enable residuals in the TPC with
0273      * TPC clusters not participating to the ACTS track fit
0274      */
0275     G4TRACKING::SC_CALIBMODE = false;
0276
0277     TpcReadoutInit(runnumber);
0278     // these lines show how to override the drift velocity and time offset values set in TpcReadoutInit
0279     // G4TPC::tpc_drift_velocity_reco = 0.0073844; // cm/ns
0280     // TpcClusterZCrossingCorrection::_vdrift = G4TPC::tpc_drift_velocity_reco;
0281     // G4TPC::tpc_tzero_reco = -5*50;  // ns
0282     std::cout << " run: " << runnumber << " samples: " << TRACKING::reco_tpc_maxtime_sample << " pre: " << TRACKING::reco_tpc_time_presample << " vdrift: " << G4TPC::tpc_drift_velocity_reco << std::endl;
0283
0284     G4TPC::REJECT_LASER_EVENTS = true;
0285     G4TPC::ENABLE_MODULE_EDGE_CORRECTIONS = true;
0286     // Flag for running the tpc hit unpacker with zero suppression on
0287     TRACKING::tpc_zero_supp = true;
0288
0289     // to turn on the default static corrections, enable the two lines below
0290     G4TPC::ENABLE_STATIC_CORRECTIONS = true;
0291     G4TPC::USE_PHI_AS_RAD_STATIC_CORRECTIONS = false;
0292
0293     // to turn on the average corrections, enable the three lines below
0294     // note: these are designed to be used only if static corrections are also applied
0295     G4TPC::ENABLE_AVERAGE_CORRECTIONS = true;
0296     G4TPC::USE_PHI_AS_RAD_AVERAGE_CORRECTIONS = false;
0297     // to use a custom file instead of the database file:
0298     G4TPC::average_correction_filename = CDBInterface::instance()->getUrl("TPC_LAMINATION_FIT_CORRECTION");
0299
0300     G4MAGNET::magfield_rescale = 1;
0301     TrackingInit();
0302
0303     output_dir = "./"; // Top dir of where the output nTuples will be written
0304     trailer = "_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + ".root";
0305
0306     if (DoUnpacking)
0307     {
0308         for (int felix = 0; felix < 6; felix++)
0309         {
0310             Mvtx_HitUnpacking(std::to_string(felix));
0311         }
0312         for (int server = 0; server < 8; server++)
0313         {
0314             Intt_HitUnpacking(std::to_string(server));
0315         }
0316         std::ostringstream ebdcname;
0317         if (runspecies == "run3pp" || runspecies == "run3auau")
0318         {
0319             for (int ebdc = 0; ebdc < 24; ebdc++)
0320             {
0321                 for (int endpoint = 0; endpoint < 2; endpoint++)
0322                 {
0323                     ebdcname.str("");
0324                     if (ebdc < 10)
0325                     {
0326                         ebdcname << "0";
0327                     }
0328                     ebdcname << ebdc << "_" << endpoint;
0329                     Tpc_HitUnpacking(ebdcname.str());
0330                 }
0331             }
0332         }
0333         if (runspecies == "run2pp")
0334         {
0335             for (int ebdc = 0; ebdc < 24; ebdc++)
0336             {
0337                 ebdcname.str("");
0338                 if (ebdc < 10)
0339                 {
0340                     ebdcname << "0";
0341                 }
0342                 ebdcname << ebdc;
0343                 Tpc_HitUnpacking(ebdcname.str());
0344             }
0345         }
0346
0347         Micromegas_HitUnpacking();
0348
0349         //se->registerSubsystem(new MvtxRawHitQA);
0350         //se->registerSubsystem(new InttRawHitQA);
0351         //se->registerSubsystem(new TpcRawHitQA);
0352     }
0353
0354     if (DoClustering)
0355     {
0356         Mvtx_Clustering();
0357
0358         Intt_Clustering();
0359
0360         Tpc_LaserEventIdentifying();
0361
0362         TPC_Clustering_run2pp();
0363
0364         Micromegas_Clustering();
0365
0366         Reject_Laser_Events();
0367
0368         //se->registerSubsystem(new MvtxClusterQA);
0369         //se->registerSubsystem(new InttClusterQA);
0370         //se->registerSubsystem(new TpcClusterQA);
0371         //se->registerSubsystem(new MicromegasClusterQA);
0372     }
0373
0374     if (DoSeeding)
0375     {
0376         Tracking_Reco_TrackSeed_run2pp();
0377
0378         auto converter = new TrackSeedTrackMapConverter("SiliconSeedConverter");
0379         // Default set to full SvtxTrackSeeds. Can be set to
0380         // SiliconTrackSeedContainer or TpcTrackSeedContainer
0381         converter->setTrackSeedName("SiliconTrackSeedContainer");
0382         converter->setTrackMapName("SiliconSvtxTrackMap");
0383         converter->setFieldMap(G4MAGNET::magfield_tracking);
0384         converter->Verbosity(0);
0385         se->registerSubsystem(converter);
0386
0387         auto finder = new PHSimpleVertexFinder("SiliconVertexFinder");
0388         finder->Verbosity(0);
0389         finder->setDcaCut(0.1);
0390         finder->setTrackPtCut(0.1);
0391         finder->setBeamLineCut(1);
0392         finder->setTrackQualityCut(1000000000);
0393         finder->setNmvtxRequired(3);
0394         finder->setOutlierPairCut(0.1);
0395         finder->setTrackMapName("SiliconSvtxTrackMap");
0396         finder->setVertexMapName("SiliconSvtxVertexMap");
0397         se->registerSubsystem(finder);
0398
0399         //auto siliconqa = new SiliconSeedsQA;
0400         //siliconqa->setTrackMapName("SiliconSvtxTrackMap");
0401         //siliconqa->setVertexMapName("SiliconSvtxVertexMap");
0402         //se->registerSubsystem(siliconqa);
0403
0404         auto convertertpc = new TrackSeedTrackMapConverter("TpcSeedConverter");
0405         // Default set to full SvtxTrackSeeds. Can be set to
0406         // SiliconTrackSeedContainer or TpcTrackSeedContainer
0407         convertertpc->setTrackSeedName("TpcTrackSeedContainer");
0408         convertertpc->setTrackMapName("TpcSvtxTrackMap");
0409         convertertpc->setFieldMap(G4MAGNET::magfield_tracking);
0410         convertertpc->Verbosity(0);
0411         se->registerSubsystem(convertertpc);
0412
0413         auto findertpc = new PHSimpleVertexFinder("TpcSimpleVertexFinder");
0414         findertpc->Verbosity(0);
0415         findertpc->setDcaCut(0.5);
0416         findertpc->setTrackPtCut(0.2);
0417         findertpc->setBeamLineCut(1);
0418         findertpc->setTrackQualityCut(1000000000);
0419         // findertpc->setNmvtxRequired(3);
0420         findertpc->setRequireMVTX(false);
0421         findertpc->setOutlierPairCut(0.1);
0422         findertpc->setTrackMapName("TpcSvtxTrackMap");
0423         findertpc->setVertexMapName("TpcSvtxVertexMap");
0424         se->registerSubsystem(findertpc);
0425
0426         //auto tpcqa = new TpcSeedsQA;
0427         //tpcqa->setTrackMapName("TpcSvtxTrackMap");
0428         //tpcqa->setVertexMapName("TpcSvtxVertexMap");
0429         //tpcqa->setSegment(rc->get_IntFlag("RUNSEGMENT"));
0430         //se->registerSubsystem(tpcqa);
0431     }
0432
0433     if (DoFitting)
0434     {
0435         if(!convertSeeds) Tracking_Reco_TrackMatching_run2pp();
0436
0437         G4TRACKING::convert_seeds_to_svtxtracks = convertSeeds;
0438         std::cout << "Converting to seeds : " << G4TRACKING::convert_seeds_to_svtxtracks << std::endl;
0439         /*
0440          * Either converts seeds to tracks with a straight line/helix fit
0441          * or run the full Acts track kalman filter fit
0442          */
0443         if (G4TRACKING::convert_seeds_to_svtxtracks)
0444         {
0445             auto *converter = new TrackSeedTrackMapConverter;
0446             // Default set to full SvtxTrackSeeds. Can be set to
0447             // SiliconTrackSeedContainer or TpcTrackSeedContainer
0448             converter->setTrackSeedName("TpcTrackSeedContainer");
0449             converter->setFieldMap(G4MAGNET::magfield_tracking);
0450             converter->useKFfit(true);
0451             converter->Verbosity(0);
0452             se->registerSubsystem(converter);
0453
0454             auto findertpc = new PHSimpleVertexFinder("TpcSimpleVertexFinder");
0455             findertpc->set_pp_mode(false);
0456             findertpc->Verbosity(0);
0457             findertpc->setDcaCut(1.5);
0458             findertpc->setTrackPtCut(0.16);
0459             findertpc->setBeamLineCut(1);
0460             findertpc->setBeamSpotCutX(-2.,2.);
0461             findertpc->setBeamSpotCutY(-2.,2.);
0462             findertpc->setTrackQualityCut(3000);
0463             // findertpc->setNmvtxRequired(3);
0464             findertpc->setRequireMVTX(false);
0465             findertpc->setOutlierPairCut(1.5);
0466             //findertpc->setTrackMapName("TpcSvtxTrackMap");
0467             //findertpc->setVertexMapName("TpcSvtxVertexMap");
0468             se->registerSubsystem(findertpc);
0469 /*
0470             auto *vtxProp = new PHActsVertexPropagator;
0471             vtxProp->Verbosity(0);
0472             vtxProp->fieldMap(G4MAGNET::magfield_tracking);
0473             se->registerSubsystem(vtxProp);
0474 */
0475         }
0476         else
0477         {
0478             Tracking_Reco_TrackFit_run2pp();
0479             //Tracking_Reco_Vertex_run2pp();
0480
0481             auto findertpc = new PHSimpleVertexFinder("TpcSimpleVertexFinder");
0482             findertpc->set_pp_mode(false);
0483             findertpc->Verbosity(2);
0484             findertpc->setDcaCut(1.5);
0485             findertpc->setTrackPtCut(0.16);
0486             findertpc->setBeamLineCut(1);
0487             findertpc->setBeamSpotCutX(-2.,2.);
0488             findertpc->setBeamSpotCutY(-2.,2.);
0489             findertpc->setTrackQualityCut(3000);
0490             // findertpc->setNmvtxRequired(3);
0491             findertpc->setRequireMVTX(false);
0492             findertpc->setOutlierPairCut(1.5);
0493             //findertpc->setTrackMapName("TpcSvtxTrackMap");
0494             //findertpc->setVertexMapName("TpcSvtxVertexMap");
0495             se->registerSubsystem(findertpc);
0496
0497             auto *vtxProp = new PHActsVertexPropagator;
0498             vtxProp->Verbosity(1);
0499             vtxProp->fieldMap(G4MAGNET::magfield_tracking);
0500             se->registerSubsystem(vtxProp);
0501         }
0502
0503         //se->registerSubsystem(new TpcSiliconQA);
0504         //se->registerSubsystem(new TrackFittingQA);
0505         //se->registerSubsystem(new VertexQA);
0506     }
0507
0508     /*
0509     if (get_calo_info)
0510     {
0511         std::cout << "Begin my calo reco" << std::endl;
0512
0513         Process_Calo_Calib();
0514
0515         Global_Reco();
0516
0517         bool doEMcalRadiusCorr = true;
0518         auto projection = new PHActsTrackProjection("CaloProjection");
0519         if (doEMcalRadiusCorr)
0520         {
0521             projection->setLayerRadius(SvtxTrack::CEMC, cemc_proj_radius);
0522         }
0523         se->registerSubsystem(projection);
0524
0525         CaloGeomMapping *cgm = new CaloGeomMapping();
0526         cgm->set_detector_name("CEMC");
0527         cgm->set_UseDetailedGeometry(true);
0528         se->registerSubsystem(cgm);
0529
0530         //////////////////
0531         // Clusters
0532         std::cout << "Building clusters" << std::endl;
0533         RawClusterBuilderTemplate *ClusterBuilder = new RawClusterBuilderTemplate("EmcRawClusterBuilderTemplate");
0534         ClusterBuilder->Detector("CEMC");
0535         ClusterBuilder->set_UseDetailedGeometry(true);
0536         ClusterBuilder->set_threshold_energy(0.070); // for when using basic calibration
0537         std::string emc_prof = getenv("CALIBRATIONROOT");
0538         emc_prof += "/EmcProfile/CEMCprof_Thresh30MeV.root";
0539         // ClusterBuilder->set_UseAltZVertex(3); //0: GlobalVertexMap, 1: MbdVertexMap, 2: Nothing, 3: G4TruthInfo
0540         ClusterBuilder->LoadProfile(emc_prof);
0541         ClusterBuilder->set_UseTowerInfo(1); // to use towerinfo objects rather than old RawTower
0542         se->registerSubsystem(ClusterBuilder);
0543     }
0544         */
0545
0546    // const std::string kappaLookup = "/sphenix/user/aopatton/MomentumScaleStudies/CalibrationsV2/kappa_lookup_79528.csv";
0547     //const std::string kappaLookup = "/sphenix/user/aopatton/MomentumScaleStudies/Calibrations_pT/kappa_lookup_79526.csv";
0548
0549     const std::string kappaLookup = "/sphenix/user/aopatton/MomentumScaleStudies/Calibrations_compare/kappa_lookup_79528_pt.csv";
0550 /*
0551     //new module here
0552     auto* momScale = new NNScaleMap("NNScaleMap");
0553
0554     momScale->setUseCDB(false);  // NNScaleMap defaults to CDB lookup; use the fixed local file instead
0555     momScale->setKappaLookupFile(kappaLookup);
0556     //momScale->setInputTrackMapName("SvtxTrackMap");
0557     //momScale->setOutputTrackMapName("SvtxTrackMapMomentumScaleCorrected");
0558     momScale->setScaleCovariance(true);
0559     se->registerSubsystem(momScale);
0560 */
0561
0562     output_dir = outDir;
0563
0564     if (run_pipi_reco)
0565         create_hf_directories(pipi_reconstruction_name, pipi_output_dir, pipi_output_reco_file);
0566     if (run_Kpi_reco)
0567         create_hf_directories(Kpi_reconstruction_name, Kpi_output_dir, Kpi_output_reco_file);
0568     if (run_KK_reco)
0569         create_hf_directories(KK_reconstruction_name, KK_output_dir, KK_output_reco_file);
0570     if (run_ppi_reco)
0571         create_hf_directories(ppi_reconstruction_name, ppi_output_dir, ppi_output_reco_file);
0572     //if (run_ee_reco)
0573     //    create_hf_directories(ee_reconstruction_name, ee_output_dir, ee_output_reco_file);
0574
0575     if (run_pipi_reco || run_Kpi_reco || run_KK_reco || run_ppi_reco)// || run_ee_reco)
0576         init_kfp_dependencies();
0577
0578     if (run_pipi_reco)
0579         reconstruct_pipi_mass();
0580     if (run_Kpi_reco)
0581         reconstruct_Kpi_mass();
0582     if (run_KK_reco)
0583         reconstruct_KK_mass();
0584     if (run_ppi_reco)
0585         reconstruct_ppi_mass();
0586     //if (run_ee_reco)
0587         //reconstruct_ee_mass();
0588
0589     se->skip(nSkip);
0590     se->run(nEvents);
0591     se->End();
0592     se->PrintTimer();
0593
0594     /*
0595     TString qaname = output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0596     std::string makeDirectory = "mkdir -p " + output_dir + "qaOut";
0597     if (run_pipi_reco)
0598     {
0599         qaname = pipi_output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0600         makeDirectory = "mkdir -p " + pipi_output_dir + "qaOut";
0601     }
0602     else if (run_Kpi_reco)
0603     {
0604         qaname = Kpi_output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0605         makeDirectory = "mkdir -p " + Kpi_output_dir + "qaOut";
0606     }
0607     else if (run_KK_reco)
0608     {
0609         qaname = KK_output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0610         makeDirectory = "mkdir -p " + KK_output_dir + "qaOut";
0611     }
0612     else if (run_ppi_reco)
0613     {
0614         qaname = ppi_output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0615         makeDirectory = "mkdir -p " + ppi_output_dir + "qaOut";
0616     }
0617     else if (run_ee_reco)
0618     {
0619         qaname = ee_output_dir + "qaOut/HIST_" + nice_runnumber.str() + "_" + nice_segment.str() + "_" + nice_skip.str() + "_QA.root";
0620         makeDirectory = "mkdir -p " + ee_output_dir + "qaOut";
0621     }
0622     std::cout << "Output QA file: " << qaname << std::endl;
0623     system(makeDirectory.c_str());
0624     std::string qaOutputFileName(qaname.Data());
0625     QA_Output(qaOutputFileName);
0626     */
0627
0628     if (run_pipi_reco)
0629         end_kfparticle(pipi_output_reco_file, pipi_output_dir);
0630     if (run_Kpi_reco)
0631         end_kfparticle(Kpi_output_reco_file, Kpi_output_dir);
0632     if (run_KK_reco)
0633         end_kfparticle(KK_output_reco_file, KK_output_dir);
0634     if (run_ppi_reco)
0635         end_kfparticle(ppi_output_reco_file, ppi_output_dir);
0636     //if (run_ee_reco)
0637     //    end_kfparticle(ee_output_reco_file, ee_output_dir);
0638
0639     delete se;
0640
0641     std::cout << "Finished" << std::endl;
0642     gSystem->Exit(0);
0643 }