ParticleGen-Jingyu/macro/Fun4All_physiTuto_Jaein.C

0001 // Fun4All headers
0002 #include <fun4all/Fun4AllServer.h>
0003 #include <fun4all/Fun4AllDstInputManager.h>
0004
0005 #include <GlobalVariables.C>
0006
0007 //#include <G4Setup_sPHENIX_Bbc.C>
0008 #include <G4Setup_sPHENIX.C>
0009 // #include <G4_Bbc.C>
0010 // #include <G4_CaloTrigger.C>
0011 // #include <G4_Centrality.C>
0012 // #include <G4_DSTReader.C>
0013 // #include <G4_Global.C>
0014 // #include <G4_HIJetReco.C>
0015 #include <G4_Input.C>
0016 // #include <G4_Jets.C>
0017 // #include <G4_KFParticle.C>
0018 // #include <G4_ParticleFlow.C>
0019 // #include <G4_Production.C>
0020 // #include <G4_TopoClusterReco.C>
0021
0022 #include <Trkr_RecoInit.C>
0023 #include <Trkr_Clustering.C>
0024 #include <Trkr_LaserClustering.C>
0025 #include <Trkr_Reco.C>
0026 #include <Trkr_Eval.C>
0027 #include <G4_ActsGeom.C>
0028 // #include <Trkr_QA.C>
0029
0030 // #include <Trkr_Diagnostics.C>
0031 // #include <G4_User.C>
0032 // #include <QA.C>
0033
0034 #include <ffamodules/FlagHandler.h>
0035 #include <ffamodules/HeadReco.h>
0036 #include <ffamodules/SyncReco.h>
0037 #include <ffamodules/CDBInterface.h>
0038
0039 #include <fun4all/Fun4AllDstOutputManager.h>
0040 #include <fun4all/Fun4AllOutputManager.h>
0041 #include <fun4all/Fun4AllServer.h>
0042
0043 #include <phool/PHRandomSeed.h>
0044 #include <phool/recoConsts.h>
0045
0046 R__LOAD_LIBRARY(libfun4all.so)
0047
0048 // #include <physiTuto/tutorial.h>
0049 #include <tutorial.h>
0050 R__LOAD_LIBRARY( libtutorial.so )
0051
0052 int Fun4All_physiTuto(
0053              int nEvents = 5,
0054              const double particle_pT = 1,
0055          const double energy_range = 0.1,
0056          const string particle_species = "PionPlus",
0057          const string &output_directory = "/sphenix/user/jaein213/INTT_Jaein/INTT/general_codes/tkumaoka/InttSeedingTrackDev/ParticleGen/output/",
0058          const string &output_filename = "results.root"
0059              )
0060 {
0061   const int skip = 0;
0062     const bool is_pythia = false;
0063   const int Nparticle = 1;
0064
0065   // const bool is_pythia = true;
0066   map<string, string> particle_map  = {
0067     {"MuonPlus", "mu+"},
0068     {"Muon", "mu-"},
0069     {"PionMinus", "pi-"},
0070     {"PionPlus", "pi+"},
0071     {"Electron", "e-"},
0072     {"Positron", "e+"},
0073     {"Proton", "proton"},
0074     {"KaonMinus", "kaon-"},
0075     {"KaonPlus", "kaon+"}
0076   };
0077
0078   if (particle_map.find(particle_species) == particle_map.end())
0079   {
0080     cout << "Invalid particle species: " << particle_species << endl;
0081     return -1;
0082   }
0083
0084   Fun4AllServer *se = Fun4AllServer::instance();
0085   //  se->Print("NODETREE"); // useless
0086   //se->Verbosity(0);
0087
0088   //Opt to print all random seed used for debugging reproducibility. Comment out to reduce stdout prints.
0089   //PHRandomSeed::Verbosity(1);
0090
0091   // just if we set some flags somewhere in this macro
0092   recoConsts *rc = recoConsts::instance();
0093   // By default every random number generator uses
0094   // PHRandomSeed() which reads /dev/urandom to get its seed
0095   // if the RANDOMSEED flag is set its value is taken as seed
0096   // You can either set this to a random value using PHRandomSeed()
0097   // which will make all seeds identical (not sure what the point of
0098   // this would be:
0099   //  rc->set_IntFlag("RANDOMSEED",PHRandomSeed());
0100   // or set it to a fixed value so you can debug your code
0101   //  rc->set_IntFlag("RANDOMSEED", 12345);
0102
0103   // INPUTREADHITS::filename[0] = inputFile;
0104
0105   // Force to use Pythia8 or GUN
0106   Input::PYTHIA8 = is_pythia;
0107   Input::GUN = !( is_pythia );
0108   Input::SIMPLE = true;
0109   if( Input::PYTHIA8 == false )
0110     {
0111       Input::GUN_NUMBER = 1; // if you need 3 of them
0112       Input::GUN_VERBOSITY = 1;
0113     }
0114
0115   // Initialize the selected Input/Event generation
0116   InputInit();
0117
0118   ///////////////
0119   // Set generator specific options
0120   ////////////////
0121   // can only be set after InputInit() is called
0122
0123   // pythia8
0124   if (Input::PYTHIA8)
0125     {
0126       //! Nominal collision geometry is selected by Input::BEAM_CONFIGURATION
0127       Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia8);
0128     }
0129   else // for GUN
0130     {
0131       // particle gun
0132       // if you run more than one of these Input::GUN_NUMBER > 1
0133       // add the settings for other with [1], next with [2]...
0134       // INPUTGENERATOR::Gun[0]->AddParticle("pi-", 0, 1, 0.1 );
0135       // INPUTGENERATOR::Gun[0]->set_vtx(0, 0, 0);
0136
0137       INPUTGENERATOR::SimpleEventGenerator[0] -> add_particles(particle_map[particle_species.c_str()], Nparticle );
0138       // INPUTGENERATOR::SimpleEventGenerator[0] -> set_name("mu-");
0139       INPUTGENERATOR::SimpleEventGenerator[0] -> set_vtx(0, 0, 0);
0140       INPUTGENERATOR::SimpleEventGenerator[0] -> set_pt_range(particle_pT - energy_range, particle_pT + energy_range);
0141       // INPUTGENERATOR::SimpleEventGenerator[0] -> set_mom(120, 0, 0);
0142       INPUTGENERATOR::SimpleEventGenerator[0] -> set_eta_range(-1, 1);
0143       // INPUTGENERATOR::SimpleEventGenerator[0] -> set_eta_range(-0.01, 0.01);
0144       INPUTGENERATOR::SimpleEventGenerator[0] -> set_phi_range(-M_PI, M_PI);
0145     }
0146
0147
0148
0149   InputRegister();
0150
0151   // Flag Handler is always needed to read flags from input (if used)
0152   // and update our rc flags with them. At the end it saves all flags
0153   // again on the DST in the Flags node under the RUN node
0154   FlagHandler *flag = new FlagHandler();
0155   se->registerSubsystem(flag);
0156
0157   //======================
0158   // Write the DST
0159   //======================
0160
0161   //Enable::DSTOUT = true;
0162   Enable::DSTOUT_COMPRESS = false;
0163   // DstOut::OutputDir = outdir;
0164   // DstOut::OutputFile = outputFile;
0165
0166   //Option to convert DST to human command readable TTree for quick poke around the outputs
0167   //  Enable::DSTREADER = true;
0168
0169   // turn the display on (default off)
0170   //Enable::DISPLAY = true;
0171
0172   //======================
0173   // What to run
0174   //======================
0175   // QA, main switch
0176   Enable::QA                    = false;
0177
0178   Enable::PIPE                  = true;
0179   Enable::PIPE_ABSORBER = true;
0180
0181   // central tracking
0182   Enable::MVTX                  = true;
0183   Enable::MVTX_CELL         = Enable::MVTX && true;
0184   Enable::MVTX_CLUSTER  = Enable::MVTX_CELL && true;
0185   Enable::MVTX_QA               = Enable::MVTX_CLUSTER && Enable::QA && true;
0186
0187   Enable::INTT                  = true;
0188   //  Enable::INTT_ABSORBER         = true; // enables layerwise support structure readout
0189   //  Enable::INTT_SUPPORT          = true; // enable global support structure readout
0190   Enable::INTT_CELL             = Enable::INTT && true;
0191   Enable::INTT_CLUSTER      = Enable::INTT_CELL && true;
0192   Enable::INTT_QA                 = Enable::INTT_CLUSTER && Enable::QA && true;
0193
0194
0195
0196   Enable::TPC                     = true;
0197   Enable::TPC_ABSORBER  = true;
0198   Enable::TPC_CELL          = Enable::TPC && true;
0199   Enable::TPC_CLUSTER       = Enable::TPC_CELL &&true;
0200   Enable::TPC_QA                = Enable::TPC_CLUSTER && Enable::QA && false;
0201
0202   Enable::MICROMEGAS                = true;
0203   Enable::MICROMEGAS_CELL           = Enable::MICROMEGAS && true;
0204   Enable::MICROMEGAS_CLUSTER    = Enable::MICROMEGAS_CELL && true;
0205   Enable::MICROMEGAS_QA             = Enable::MICROMEGAS_CLUSTER && Enable::QA && false;
0206
0207   Enable::TRACKING_TRACK    =  true;
0208
0209   Enable::TRACKING_EVAL     = Enable::TRACKING_TRACK && false;
0210   Enable::TRACKING_QA           = Enable::TRACKING_TRACK && Enable::QA && true; // <-- this one
0211
0212
0213   //! forward flux return plug door. Out of acceptance and off by default.
0214   //Enable::PLUGDOOR                = true;
0215   Enable::PLUGDOOR_ABSORBER         = true;
0216
0217
0218   Enable::CEMC = true;
0219   Enable::CEMC_ABSORBER = true;
0220   Enable::CEMC_CELL = Enable::CEMC && true;
0221   Enable::CEMC_TOWER = Enable::CEMC_CELL && true;
0222   Enable::CEMC_CLUSTER = Enable::CEMC_TOWER && true;
0223   Enable::CEMC_EVAL = Enable::CEMC_CLUSTER && true;
0224   Enable::CEMC_QA = Enable::CEMC_CLUSTER && Enable::QA && false;
0225
0226   Enable::HCALIN = true;
0227   Enable::HCALIN_ABSORBER = true;
0228   Enable::HCALIN_CELL = Enable::HCALIN && true;
0229   Enable::HCALIN_TOWER = Enable::HCALIN_CELL && true;
0230   Enable::HCALIN_CLUSTER = Enable::HCALIN_TOWER && true;
0231   Enable::HCALIN_EVAL = Enable::HCALIN_CLUSTER && true;
0232   Enable::HCALIN_QA = Enable::HCALIN_CLUSTER && Enable::QA && false;
0233
0234
0235   Enable::MAGNET                = true;
0236   Enable::MAGNET_ABSORBER           = true;
0237
0238
0239   Enable::HCALOUT = true;
0240   Enable::HCALOUT_ABSORBER = true;
0241   Enable::HCALOUT_CELL = Enable::HCALOUT && true;
0242   Enable::HCALOUT_TOWER = Enable::HCALOUT_CELL && true;
0243   Enable::HCALOUT_CLUSTER = Enable::HCALOUT_TOWER && true;
0244   Enable::HCALOUT_EVAL = Enable::HCALOUT_CLUSTER && true;
0245   Enable::HCALOUT_QA = Enable::HCALOUT_CLUSTER && Enable::QA && false;
0246
0247
0248   // new settings using Enable namespace in GlobalVariables.C
0249   Enable::BLACKHOLE             = true;
0250   //Enable::BLACKHOLE_SAVEHITS          = false; // turn off saving of bh hits
0251   //Enable::BLACKHOLE_FORWARD_SAVEHITS      = false; // disable forward/backward hits
0252   //BlackHoleGeometry::visible          = true;
0253
0254   //////////////////////////////
0255   // conditions DB flags      //
0256   //////////////////////////////
0257   Enable::CDB                   = true;
0258   // global tag
0259   rc->set_StringFlag("CDB_GLOBALTAG",CDB::global_tag);
0260   // 64 bit timestamp
0261   rc->set_uint64Flag("TIMESTAMP",CDB::timestamp);
0262
0263   /////////////////
0264   // Magnet Settings
0265   /////////////////
0266   G4MAGNET::magfield_rescale = 1.;  // for zero field
0267
0268   // Initialize the selected subsystems
0269   G4Init();
0270
0271   ///////////////////////
0272   // GEANT4 Detector description
0273   ///////////////////////
0274   if (!Input::READHITS) G4Setup();
0275
0276   ///////////////////////
0277   // Detector Division //
0278   ///////////////////////
0279
0280   //  if ((Enable::MBD && Enable::MBDRECO) || Enable::MBDFAKE) Mbd_Reco();
0281   if (Enable::MVTX_CELL) Mvtx_Cells();
0282   if (Enable::INTT_CELL) Intt_Cells();
0283
0284   //////////////////////////////////
0285   // CEMC towering and clustering //
0286   //////////////////////////////////
0287   if (Enable::CEMC_CELL) CEMC_Cells();
0288   if (Enable::CEMC_TOWER) CEMC_Towers();
0289   if (Enable::CEMC_CLUSTER) CEMC_Clusters();
0290
0291   if (Enable::HCALIN_CELL) HCALInner_Cells();
0292   if (Enable::HCALOUT_CELL) HCALOuter_Cells();
0293   if (Enable::HCALIN_TOWER) HCALInner_Towers();
0294   if (Enable::HCALIN_CLUSTER) HCALInner_Clusters();
0295   if (Enable::HCALOUT_TOWER) HCALOuter_Towers();
0296   if (Enable::HCALOUT_CLUSTER) HCALOuter_Clusters();
0297
0298
0299   /////////////////
0300   // SVTX tracking
0301   //////////////--
0302   if(Enable::TRACKING_TRACK) TrackingInit();
0303   if (Enable::MVTX_CLUSTER) Mvtx_Clustering();
0304   if (Enable::INTT_CLUSTER) Intt_Clustering();
0305   if (Enable::TRACKING_TRACK) Tracking_Reco();
0306
0307   //////////////////////
0308   // Simulation evaluation
0309   //////////////////////
0310   // string outputroot = outputFile;
0311   // string remove_this = ".root";
0312   // size_t pos = outputroot.find(remove_this);
0313   // if (pos != string::npos) outputroot.erase(pos, remove_this.length());
0314
0315   // if (Enable::TRACKING_EVAL) Tracking_Eval(outputroot + "_g4svtx_eval.root");
0316
0317   tutorial* analysis_module = new tutorial(
0318     "name",
0319     output_directory,
0320     output_filename
0321   );
0322   // string output_path = "tutorial_results_";
0323   // output_path += "pythia8_MC.root";
0324
0325   // analysis_module->SetOutputPath( output_path );
0326   se->registerSubsystem( analysis_module );
0327
0328   // Flag Handler is always needed to read flags from input (if used)
0329   // and update our rc flags with them. At the end it saves all flags
0330   // again on the DST in the Flags node under the RUN node
0331   // FlagHandler *flag = new FlagHandler();
0332   // se->registerSubsystem(flag);
0333   se->skip(skip);
0334   se->run(nEvents);
0335   se->End();
0336   delete se;
0337
0338   gSystem->Exit(0);
0339   return 0;
0340 }
0341
0342