macro/ana.354/Fun4All_G4_sPHENIX.C

0001 #ifndef MACRO_FUN4ALLG4SPHENIX_C
0002 #define MACRO_FUN4ALLG4SPHENIX_C
0003
0004 // c++ includes --
0005 #include <string>
0006 #include <iostream>
0007 #include <fstream>
0008 #include <math.h>
0009 #include <set>
0010 #include <map>
0011
0012 using std::cout;
0013 using std::endl;
0014 using std::string;
0015 using std::to_string;
0016 using std::set;
0017 using std::istringstream;
0018 using std::max;
0019 using std::map;
0020
0021 // root includes --
0022 #include <TSystem.h>
0023 #include <TROOT.h>
0024
0025 #include <GlobalVariables.C>
0026
0027 #include <DisplayOn.C>
0028 #include "G4Setup_sPHENIX.C"
0029 #include <G4_Bbc.C>
0030 #include <G4_CaloTrigger.C>
0031 #include <G4_Centrality.C>
0032 #include <G4_DSTReader.C>
0033 #include <G4_Global.C>
0034 #include <G4_HIJetReco.C>
0035 #include <G4_Input.C>
0036 #include <G4_Jets.C>
0037 #include <G4_KFParticle.C>
0038 #include <G4_ParticleFlow.C>
0039 #include <G4_Production.C>
0040 #include <G4_TopoClusterReco.C>
0041 #include <G4_Tracking.C>
0042 #include <G4_TrackingDiagnostics.C>
0043 #include <G4_User.C>
0044 #include <QA.C>
0045
0046 #include <ffamodules/FlagHandler.h>
0047 #include <ffamodules/HeadReco.h>
0048 #include <ffamodules/SyncReco.h>
0049 #include <ffamodules/XploadInterface.h>
0050
0051 #include <fun4all/Fun4AllDstOutputManager.h>
0052 #include <fun4all/Fun4AllOutputManager.h>
0053 #include <fun4all/Fun4AllServer.h>
0054
0055 #include <phool/PHRandomSeed.h>
0056 #include <phool/recoConsts.h>
0057
0058 R__LOAD_LIBRARY(libfun4all.so)
0059 R__LOAD_LIBRARY(libffamodules.so)
0060
0061 // For HepMC Hijing
0062 // try inputFile = /sphenix/sim/sim01/sphnxpro/sHijing_HepMC/sHijing_0-12fm.dat
0063
0064 int Fun4All_G4_sPHENIX(
0065     const int nEvents = 1,
0066     const int seed = 0,
0067     const string &outdir = ".",
0068     const string &tag = "0",
0069     const string &outputFile = "G4sPHENIX.root",
0070     const string &inputFile = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root",
0071     const string &embed_input_file = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root",
0072     const int skip = 0)
0073 {
0074   Fun4AllServer *se = Fun4AllServer::instance();
0075   se->Verbosity(1);
0076
0077   //Opt to print all random seed used for debugging reproducibility. Comment out to reduce stdout prints.
0078   PHRandomSeed::Verbosity(1);
0079
0080   // just if we set some flags somewhere in this macro
0081   recoConsts *rc = recoConsts::instance();
0082   // By default every random number generator uses
0083   // PHRandomSeed() which reads /dev/urandom to get its seed
0084   // if the RANDOMSEED flag is set its value is taken as seed
0085   // You can either set this to a random value using PHRandomSeed()
0086   // which will make all seeds identical (not sure what the point of
0087   // this would be:
0088   //  rc->set_IntFlag("RANDOMSEED",PHRandomSeed());
0089   // or set it to a fixed value so you can debug your code
0090   //  rc->set_IntFlag("RANDOMSEED", 12345);
0091
0092   if(seed > 0) {
0093     rc->set_IntFlag("RANDOMSEED", seed);
0094   }
0095
0096   //===============
0097   // Input options
0098   //===============
0099   // verbosity setting (applies to all input managers)
0100   Input::VERBOSITY = 0;
0101   // First enable the input generators
0102   // Either:
0103   // read previously generated g4-hits files, in this case it opens a DST and skips
0104   // the simulations step completely. The G4Setup macro is only loaded to get information
0105   // about the number of layers used for the cell reco code
0106   //  Input::READHITS = true;
0107   INPUTREADHITS::filename[0] = inputFile;
0108   // if you use a filelist
0109   // INPUTREADHITS::listfile[0] = inputFile;
0110   // Or:
0111   // Use particle generator
0112   // And
0113   // Further choose to embed newly simulated events to a previous simulation. Not compatible with `readhits = true`
0114   // In case embedding into a production output, please double check your G4Setup_sPHENIX.C and G4_*.C consistent with those in the production macro folder
0115   // E.g. /sphenix/sim//sim01/production/2016-07-21/single_particle/spacal2d/
0116   //  Input::EMBED = true;
0117   INPUTEMBED::filename[0] = embed_input_file;
0118   // if you use a filelist
0119   //INPUTEMBED::listfile[0] = embed_input_file;
0120
0121   Input::SIMPLE = true;
0122   // Input::SIMPLE_NUMBER = 2; // if you need 2 of them
0123   // Input::SIMPLE_VERBOSITY = 1;
0124
0125   // Enable this is emulating the nominal pp/pA/AA collision vertex distribution
0126   // Input::BEAM_CONFIGURATION = Input::AA_COLLISION; // Input::AA_COLLISION (default), Input::pA_COLLISION, Input::pp_COLLISION
0127
0128   //  Input::PYTHIA6 = true;
0129
0130   // Input::PYTHIA8 = true;
0131
0132   //  Input::GUN = true;
0133   //  Input::GUN_NUMBER = 3; // if you need 3 of them
0134   // Input::GUN_VERBOSITY = 1;
0135
0136   //D0 generator
0137   //Input::DZERO = false;
0138   //Input::DZERO_VERBOSITY = 0;
0139   //Lambda_c generator //Not ready yet
0140   //Input::LAMBDAC = false;
0141   //Input::LAMBDAC_VERBOSITY = 0;
0142   // Upsilon generator
0143   //Input::UPSILON = true;
0144   //Input::UPSILON_NUMBER = 3; // if you need 3 of them
0145   //Input::UPSILON_VERBOSITY = 0;
0146
0147   //  Input::HEPMC = true;
0148   INPUTHEPMC::filename = inputFile;
0149
0150   // Event pile up simulation with collision rate in Hz MB collisions.
0151   //Input::PILEUPRATE = 100e3;
0152
0153   //-----------------
0154   // Initialize the selected Input/Event generation
0155   //-----------------
0156   // This creates the input generator(s)
0157   InputInit();
0158
0159   //--------------
0160   // Set generator specific options
0161   //--------------
0162   // can only be set after InputInit() is called
0163
0164   // Simple Input generator:
0165   // if you run more than one of these Input::SIMPLE_NUMBER > 1
0166   // add the settings for other with [1], next with [2]...
0167   if (Input::SIMPLE)
0168   {
0169     // particle id: https://pdg.lbl.gov/2020/reviews/rpp2020-rev-monte-carlo-numbering.pdf
0170     // photon id: 22
0171     INPUTGENERATOR::SimpleEventGenerator[0]->add_particles(22, 1);
0172     if (Input::HEPMC || Input::EMBED)
0173     {
0174       INPUTGENERATOR::SimpleEventGenerator[0]->set_reuse_existing_vertex(true);
0175       INPUTGENERATOR::SimpleEventGenerator[0]->set_existing_vertex_offset_vector(0.0, 0.0, 0.0);
0176     }
0177     else
0178     {
0179       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_function(PHG4SimpleEventGenerator::Gaus,
0180                                                                                 PHG4SimpleEventGenerator::Gaus,
0181                                                                                 PHG4SimpleEventGenerator::Gaus);
0182       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_mean(0., 0., 0.);
0183       INPUTGENERATOR::SimpleEventGenerator[0]->set_vertex_distribution_width(0.01, 0.01, 5.);
0184     }
0185     // generation in one sector
0186     INPUTGENERATOR::SimpleEventGenerator[0]->set_eta_range(0, 1.152);
0187     INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(-M_PI/32.0 - 0.02, M_PI/32.0 + 0.02);
0188
0189     // generation along the bottom boundary of the sector
0190     // INPUTGENERATOR::SimpleEventGenerator[0]->set_eta_range(0, 1.152);
0191     // INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(-M_PI/32.0 - 0.02, -0.073 + 0.02);
0192
0193     // generation along the top boundary of the sector
0194     // INPUTGENERATOR::SimpleEventGenerator[0]->set_eta_range(0, 1.152);
0195     // INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(0.073 - 0.02, M_PI/32.0 + 0.02);
0196
0197     // for testing use full range, all sectors
0198     // INPUTGENERATOR::SimpleEventGenerator[0]->set_eta_range(-1.152, 1.152);
0199     // INPUTGENERATOR::SimpleEventGenerator[0]->set_phi_range(-M_PI, M_PI);
0200
0201     // Train
0202     INPUTGENERATOR::SimpleEventGenerator[0]->set_pt_range(20, 21.);
0203
0204     // Test
0205     // INPUTGENERATOR::SimpleEventGenerator[0]->set_pt_range(1, 21.);
0206   }
0207   // Upsilons
0208   // if you run more than one of these Input::UPSILON_NUMBER > 1
0209   // add the settings for other with [1], next with [2]...
0210   if (Input::UPSILON)
0211   {
0212     INPUTGENERATOR::VectorMesonGenerator[0]->add_decay_particles("e", 0);
0213     INPUTGENERATOR::VectorMesonGenerator[0]->set_rapidity_range(-1, 1);
0214     INPUTGENERATOR::VectorMesonGenerator[0]->set_pt_range(0., 10.);
0215     // Y species - select only one, last one wins
0216     INPUTGENERATOR::VectorMesonGenerator[0]->set_upsilon_1s();
0217     if (Input::HEPMC || Input::EMBED)
0218     {
0219       INPUTGENERATOR::VectorMesonGenerator[0]->set_reuse_existing_vertex(true);
0220       INPUTGENERATOR::VectorMesonGenerator[0]->set_existing_vertex_offset_vector(0.0, 0.0, 0.0);
0221     }
0222   }
0223   // particle gun
0224   // if you run more than one of these Input::GUN_NUMBER > 1
0225   // add the settings for other with [1], next with [2]...
0226   if (Input::GUN)
0227   {
0228     INPUTGENERATOR::Gun[0]->AddParticle("pi-", 0, 1, 0);
0229     INPUTGENERATOR::Gun[0]->set_vtx(0, 0, 0);
0230   }
0231
0232   // pythia6
0233   if (Input::PYTHIA6)
0234   {
0235     //! Nominal collision geometry is selected by Input::BEAM_CONFIGURATION
0236     Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia6);
0237   }
0238   // pythia8
0239   if (Input::PYTHIA8)
0240   {
0241     //! Nominal collision geometry is selected by Input::BEAM_CONFIGURATION
0242     Input::ApplysPHENIXBeamParameter(INPUTGENERATOR::Pythia8);
0243   }
0244
0245   //--------------
0246   // Set Input Manager specific options
0247   //--------------
0248   // can only be set after InputInit() is called
0249
0250   if (Input::HEPMC)
0251   {
0252     //! Nominal collision geometry is selected by Input::BEAM_CONFIGURATION
0253     Input::ApplysPHENIXBeamParameter(INPUTMANAGER::HepMCInputManager);
0254
0255     // optional overriding beam parameters
0256     //INPUTMANAGER::HepMCInputManager->set_vertex_distribution_width(100e-4, 100e-4, 8, 0);  //optional collision smear in space, time
0257     //    INPUTMANAGER::HepMCInputManager->set_vertex_distribution_mean(0,0,0,0);//optional collision central position shift in space, time
0258     // //optional choice of vertex distribution function in space, time
0259     //INPUTMANAGER::HepMCInputManager->set_vertex_distribution_function(PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus, PHHepMCGenHelper::Gaus);
0260     //! embedding ID for the event
0261     //! positive ID is the embedded event of interest, e.g. jetty event from pythia
0262     //! negative IDs are backgrounds, .e.g out of time pile up collisions
0263     //! Usually, ID = 0 means the primary Au+Au collision background
0264     //INPUTMANAGER::HepMCInputManager->set_embedding_id(Input::EmbedID);
0265     if (Input::PILEUPRATE > 0)
0266     {
0267       // Copy vertex settings from foreground hepmc input
0268       INPUTMANAGER::HepMCPileupInputManager->CopyHelperSettings(INPUTMANAGER::HepMCInputManager);
0269       // and then modify the ones you want to be different
0270       // INPUTMANAGER::HepMCPileupInputManager->set_vertex_distribution_width(100e-4,100e-4,8,0);
0271     }
0272   }
0273   if (Input::PILEUPRATE > 0)
0274   {
0275     //! Nominal collision geometry is selected by Input::BEAM_CONFIGURATION
0276     Input::ApplysPHENIXBeamParameter(INPUTMANAGER::HepMCPileupInputManager);
0277   }
0278   // register all input generators with Fun4All
0279   InputRegister();
0280
0281   if (! Input::READHITS)
0282   {
0283     rc->set_IntFlag("RUNNUMBER",1);
0284
0285     SyncReco *sync = new SyncReco();
0286     se->registerSubsystem(sync);
0287
0288     HeadReco *head = new HeadReco();
0289     se->registerSubsystem(head);
0290   }
0291 // Flag Handler is always needed to read flags from input (if used)
0292 // and update our rc flags with them. At the end it saves all flags
0293 // again on the DST in the Flags node under the RUN node
0294   FlagHandler *flag = new FlagHandler();
0295   se->registerSubsystem(flag);
0296
0297   // set up production relatedstuff
0298   //   Enable::PRODUCTION = true;
0299
0300   //======================
0301   // Write the DST
0302   //======================
0303
0304   //Enable::DSTOUT = true;
0305   Enable::DSTOUT_COMPRESS = false;
0306   DstOut::OutputDir = outdir;
0307   DstOut::OutputFile = outputFile;
0308
0309   //Option to convert DST to human command readable TTree for quick poke around the outputs
0310   //  Enable::DSTREADER = true;
0311
0312   // turn the display on (default off)
0313    //Enable::DISPLAY = true;
0314
0315   //======================
0316   // What to run
0317   //======================
0318
0319   // QA, main switch
0320   Enable::QA = false;
0321
0322   // Global options (enabled for all enables subsystems - if implemented)
0323   //  Enable::ABSORBER = true;
0324   //  Enable::OVERLAPCHECK = true;
0325   //  Enable::VERBOSITY = 1;
0326
0327   // Enable::BBC = true;
0328   // Enable::BBC_SUPPORT = true; // save hist in bbc support structure
0329   Enable::BBCFAKE = true;  // Smeared vtx and t0, use if you don't want real BBC in simulation
0330
0331   Enable::PIPE = true;
0332   Enable::PIPE_ABSORBER = true;
0333
0334   // central tracking
0335   Enable::MVTX = false;
0336   Enable::MVTX_CELL = Enable::MVTX && true;
0337   Enable::MVTX_CLUSTER = Enable::MVTX_CELL && true;
0338   Enable::MVTX_QA = Enable::MVTX_CLUSTER && Enable::QA && true;
0339
0340   Enable::INTT = false;
0341 //  Enable::INTT_ABSORBER = true; // enables layerwise support structure readout
0342 //  Enable::INTT_SUPPORT = true; // enable global support structure readout
0343   Enable::INTT_CELL = Enable::INTT && true;
0344   Enable::INTT_CLUSTER = Enable::INTT_CELL && true;
0345   Enable::INTT_QA = Enable::INTT_CLUSTER && Enable::QA && true;
0346
0347   Enable::TPC = false;
0348   Enable::TPC_ABSORBER = true;
0349   Enable::TPC_CELL = Enable::TPC && true;
0350   Enable::TPC_CLUSTER = Enable::TPC_CELL && true;
0351   Enable::TPC_QA = Enable::TPC_CLUSTER && Enable::QA && true;
0352
0353   Enable::MICROMEGAS = false;
0354   Enable::MICROMEGAS_CELL = Enable::MICROMEGAS && true;
0355   Enable::MICROMEGAS_CLUSTER = Enable::MICROMEGAS_CELL && true;
0356   Enable::MICROMEGAS_QA = Enable::MICROMEGAS_CLUSTER && Enable::QA && true;
0357
0358   Enable::TRACKING_TRACK = (Enable::MICROMEGAS_CLUSTER && Enable::TPC_CLUSTER && Enable::INTT_CLUSTER && Enable::MVTX_CLUSTER) && true;
0359   Enable::TRACKING_EVAL = Enable::TRACKING_TRACK && true;
0360   Enable::TRACKING_QA = Enable::TRACKING_TRACK && Enable::QA && true;
0361
0362   //Additional tracking tools
0363   //Enable::TRACKING_DIAGNOSTICS = Enable::TRACKING_TRACK && true;
0364   //G4TRACKING::filter_conversion_electrons = true;
0365
0366
0367   //  cemc electronics + thin layer of W-epoxy to get albedo from cemc
0368   //  into the tracking, cannot run together with CEMC
0369   //  Enable::CEMCALBEDO = true;
0370
0371   Enable::CEMC = true;
0372   Enable::CEMC_ABSORBER = true;
0373   Enable::CEMC_CELL = Enable::CEMC && true;
0374   Enable::CEMC_TOWER = Enable::CEMC_CELL && true;
0375   Enable::CEMC_CLUSTER = Enable::CEMC_TOWER && true;
0376   Enable::CEMC_EVAL = Enable::CEMC_CLUSTER && true;
0377   Enable::CEMC_QA = Enable::CEMC_CLUSTER && Enable::QA && true;
0378
0379   Enable::HCALIN = false;
0380   Enable::HCALIN_ABSORBER = true;
0381   Enable::HCALIN_CELL = Enable::HCALIN && true;
0382   Enable::HCALIN_TOWER = Enable::HCALIN_CELL && true;
0383   Enable::HCALIN_CLUSTER = Enable::HCALIN_TOWER && true;
0384   Enable::HCALIN_EVAL = Enable::HCALIN_CLUSTER && true;
0385   Enable::HCALIN_QA = Enable::HCALIN_CLUSTER && Enable::QA && true;
0386
0387   Enable::MAGNET = true;
0388   Enable::MAGNET_ABSORBER = true;
0389
0390   Enable::HCALOUT = false;
0391   Enable::HCALOUT_ABSORBER = true;
0392   Enable::HCALOUT_CELL = Enable::HCALOUT && true;
0393   Enable::HCALOUT_TOWER = Enable::HCALOUT_CELL && true;
0394   Enable::HCALOUT_CLUSTER = Enable::HCALOUT_TOWER && true;
0395   Enable::HCALOUT_EVAL = Enable::HCALOUT_CLUSTER && true;
0396   Enable::HCALOUT_QA = Enable::HCALOUT_CLUSTER && Enable::QA && true;
0397
0398   Enable::EPD = false;
0399   Enable::EPD_TILE = Enable::EPD && true;
0400
0401   Enable::BEAMLINE = true;
0402 //  Enable::BEAMLINE_ABSORBER = true;  // makes the beam line magnets sensitive volumes
0403 //  Enable::BEAMLINE_BLACKHOLE = true; // turns the beamline magnets into black holes
0404   Enable::ZDC = true;
0405 //  Enable::ZDC_ABSORBER = true;
0406 //  Enable::ZDC_SUPPORT = true;
0407   Enable::ZDC_TOWER = Enable::ZDC && true;
0408   Enable::ZDC_EVAL = Enable::ZDC_TOWER && true;
0409
0410   //! forward flux return plug door. Out of acceptance and off by default.
0411   //Enable::PLUGDOOR = true;
0412   Enable::PLUGDOOR_ABSORBER = true;
0413
0414   Enable::GLOBAL_RECO = (Enable::BBCFAKE || Enable::TRACKING_TRACK) && true;
0415   //Enable::GLOBAL_FASTSIM = true;
0416
0417   //Enable::KFPARTICLE = true;
0418   //Enable::KFPARTICLE_VERBOSITY = 1;
0419   //Enable::KFPARTICLE_TRUTH_MATCH = true;
0420   //Enable::KFPARTICLE_SAVE_NTUPLE = true;
0421
0422   Enable::CALOTRIGGER = Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER && false;
0423
0424   Enable::JETS = (Enable::GLOBAL_RECO || Enable::GLOBAL_FASTSIM) && false;
0425   Enable::JETS_EVAL = Enable::JETS && true;
0426   Enable::JETS_QA = Enable::JETS && Enable::QA && true;
0427
0428   // HI Jet Reco for p+Au / Au+Au collisions (default is false for
0429   // single particle / p+p-only simulations, or for p+Au / Au+Au
0430   // simulations which don't particularly care about jets)
0431   Enable::HIJETS = Enable::JETS && Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER && false;
0432
0433   // 3-D topoCluster reconstruction, potentially in all calorimeter layers
0434   Enable::TOPOCLUSTER = Enable::CEMC_TOWER && Enable::HCALIN_TOWER && Enable::HCALOUT_TOWER && false;
0435   // particle flow jet reconstruction - needs topoClusters!
0436   Enable::PARTICLEFLOW = Enable::TOPOCLUSTER && true;
0437   // centrality reconstruction
0438   Enable::CENTRALITY = false;
0439
0440   // new settings using Enable namespace in GlobalVariables.C
0441   Enable::BLACKHOLE = false;
0442   //Enable::BLACKHOLE_SAVEHITS = false; // turn off saving of bh hits
0443   //Enable::BLACKHOLE_FORWARD_SAVEHITS = false; // disable forward/backward hits
0444   //BlackHoleGeometry::visible = true;
0445
0446   // run user provided code (from local G4_User.C)
0447   //Enable::USER = true;
0448
0449   //===============
0450   // conditions DB flags
0451   //===============
0452   //Enable::XPLOAD = true;
0453   // tag
0454   rc->set_StringFlag("XPLOAD_TAG",XPLOAD::tag);
0455   // database config
0456   rc->set_StringFlag("XPLOAD_CONFIG",XPLOAD::config);
0457   // 64 bit timestamp
0458   rc->set_uint64Flag("TIMESTAMP",XPLOAD::timestamp);
0459
0460   //---------------
0461   // World Settings
0462   //---------------
0463   //  G4WORLD::PhysicsList = "FTFP_BERT"; //FTFP_BERT_HP best for calo
0464   //  G4WORLD::WorldMaterial = "G4_AIR"; // set to G4_GALACTIC for material scans
0465
0466   //---------------
0467   // Magnet Settings
0468   //---------------
0469
0470   //  G4MAGNET::magfield =  string(getenv("CALIBRATIONROOT"))+ string("/Field/Map/sphenix3dbigmapxyz.root");  // default map from the calibration database
0471   //  G4MAGNET::magfield = "1.5"; // alternatively to specify a constant magnetic field, give a float number, which will be translated to solenoidal field in T, if string use as fieldmap name (including path)
0472 //  G4MAGNET::magfield_rescale = 1.;  // make consistent with expected Babar field strength of 1.4T
0473
0474   //---------------
0475   // Pythia Decayer
0476   //---------------
0477   // list of decay types in
0478   // $OFFLINE_MAIN/include/g4decayer/EDecayType.hh
0479   // default is All:
0480   // G4P6DECAYER::decayType = EDecayType::kAll;
0481
0482   // Initialize the selected subsystems
0483   G4Init();
0484
0485   //---------------------
0486   // GEANT4 Detector description
0487   //---------------------
0488   if (!Input::READHITS)
0489   {
0490     G4Setup();
0491   }
0492
0493   //------------------
0494   // Detector Division
0495   //------------------
0496
0497   if (Enable::BBC || Enable::BBCFAKE) Bbc_Reco();
0498
0499   if (Enable::MVTX_CELL) Mvtx_Cells();
0500   if (Enable::INTT_CELL) Intt_Cells();
0501   if (Enable::TPC_CELL) TPC_Cells();
0502   if (Enable::MICROMEGAS_CELL) Micromegas_Cells();
0503
0504   if (Enable::CEMC_CELL) CEMC_Cells();
0505
0506   if (Enable::HCALIN_CELL) HCALInner_Cells();
0507
0508   if (Enable::HCALOUT_CELL) HCALOuter_Cells();
0509
0510   //-----------------------------
0511   // CEMC towering and clustering
0512   //-----------------------------
0513
0514   if (Enable::CEMC_TOWER) CEMC_Towers();
0515   if (Enable::CEMC_CLUSTER) CEMC_Clusters();
0516
0517   //--------------
0518   // EPD tile reconstruction
0519   //--------------
0520
0521   if (Enable::EPD_TILE) EPD_Tiles();
0522
0523   //-----------------------------
0524   // HCAL towering and clustering
0525   //-----------------------------
0526
0527   if (Enable::HCALIN_TOWER) HCALInner_Towers();
0528   if (Enable::HCALIN_CLUSTER) HCALInner_Clusters();
0529
0530   if (Enable::HCALOUT_TOWER) HCALOuter_Towers();
0531   if (Enable::HCALOUT_CLUSTER) HCALOuter_Clusters();
0532
0533   // if enabled, do topoClustering early, upstream of any possible jet reconstruction
0534   if (Enable::TOPOCLUSTER) TopoClusterReco();
0535
0536   //--------------
0537   // SVTX tracking
0538   //--------------
0539   if(Enable::TRACKING_TRACK)
0540     {
0541       TrackingInit();
0542     }
0543   if (Enable::MVTX_CLUSTER) Mvtx_Clustering();
0544   if (Enable::INTT_CLUSTER) Intt_Clustering();
0545   if (Enable::TPC_CLUSTER) TPC_Clustering();
0546   if (Enable::MICROMEGAS_CLUSTER) Micromegas_Clustering();
0547
0548   if (Enable::TRACKING_TRACK)
0549   {
0550     Tracking_Reco();
0551   }
0552
0553   if(Enable::TRACKING_DIAGNOSTICS)
0554     {
0555       const std::string kshortFile = "./kshort_" + outputFile;
0556       const std::string residualsFile = "./residuals_" + outputFile;
0557
0558       G4KshortReconstruction(kshortFile);
0559       seedResiduals(residualsFile);
0560     }
0561
0562   //-----------------
0563   // Global Vertexing
0564   //-----------------
0565
0566   if (Enable::GLOBAL_RECO && Enable::GLOBAL_FASTSIM)
0567   {
0568     cout << "You can only enable Enable::GLOBAL_RECO or Enable::GLOBAL_FASTSIM, not both" << endl;
0569     gSystem->Exit(1);
0570   }
0571   if (Enable::GLOBAL_RECO)
0572   {
0573     Global_Reco();
0574   }
0575   else if (Enable::GLOBAL_FASTSIM)
0576   {
0577     Global_FastSim();
0578   }
0579
0580   //-----------------
0581   // Centrality Determination
0582   //-----------------
0583
0584   if (Enable::CENTRALITY)
0585   {
0586       Centrality();
0587   }
0588
0589   //-----------------
0590   // Calo Trigger Simulation
0591   //-----------------
0592
0593   if (Enable::CALOTRIGGER)
0594   {
0595     CaloTrigger_Sim();
0596   }
0597
0598   //---------
0599   // Jet reco
0600   //---------
0601
0602   if (Enable::JETS) Jet_Reco();
0603   if (Enable::HIJETS) HIJetReco();
0604
0605   if (Enable::PARTICLEFLOW) ParticleFlow();
0606
0607   //----------------------
0608   // Simulation evaluation
0609   //----------------------
0610   string outputroot = outdir + "/" + outputFile;
0611   string remove_this = ".root";
0612   size_t pos = outputroot.find(remove_this);
0613   if (pos != string::npos)
0614   {
0615     outputroot.erase(pos, remove_this.length());
0616   }
0617
0618   if (Enable::TRACKING_EVAL) Tracking_Eval(outputroot + "_g4svtx_eval.root");
0619
0620   if (Enable::CEMC_EVAL) CEMC_Eval(outputroot + "_g4cemc_eval-" + tag + ".root");
0621
0622   if (Enable::HCALIN_EVAL) HCALInner_Eval(outputroot + "_g4hcalin_eval.root");
0623
0624   if (Enable::HCALOUT_EVAL) HCALOuter_Eval(outputroot + "_g4hcalout_eval.root");
0625
0626   if (Enable::JETS_EVAL) Jet_Eval(outputroot + "_g4jet_eval.root");
0627
0628   if (Enable::DSTREADER) G4DSTreader(outputroot + "_DSTReader.root");
0629
0630   if (Enable::USER) UserAnalysisInit();
0631
0632   // Writes electrons from conversions to a new track map on the node tree
0633   // the ntuple file is for diagnostics, it is produced only if the flag is set in G4_Tracking.C
0634   if(G4TRACKING::filter_conversion_electrons) Filter_Conversion_Electrons(outputroot + "_secvert_ntuple.root");
0635
0636   //======================
0637   // Run KFParticle on evt
0638   //======================
0639   if (Enable::KFPARTICLE && Input::UPSILON) KFParticle_Upsilon_Reco();
0640   if (Enable::KFPARTICLE && Input::DZERO) KFParticle_D0_Reco();
0641
0642   //----------------------
0643   // Standard QAs
0644   //----------------------
0645
0646   if (Enable::CEMC_QA) CEMC_QA();
0647   if (Enable::HCALIN_QA) HCALInner_QA();
0648   if (Enable::HCALOUT_QA) HCALOuter_QA();
0649
0650   if (Enable::JETS_QA) Jet_QA();
0651
0652   if (Enable::MVTX_QA) Mvtx_QA();
0653   if (Enable::INTT_QA) Intt_QA();
0654   if (Enable::TPC_QA) TPC_QA();
0655   if (Enable::MICROMEGAS_QA) Micromegas_QA();
0656   if (Enable::TRACKING_QA) Tracking_QA();
0657
0658   if (Enable::TRACKING_QA && Enable::CEMC_QA && Enable::HCALIN_QA && Enable::HCALOUT_QA) QA_G4CaloTracking();
0659
0660   //--------------
0661   // Set up Input Managers
0662   //--------------
0663
0664   InputManagers();
0665
0666   if (Enable::PRODUCTION)
0667   {
0668     Production_CreateOutputDir();
0669   }
0670
0671   if (Enable::DSTOUT)
0672   {
0673     string FullOutFile = DstOut::OutputDir + "/" + DstOut::OutputFile;
0674     Fun4AllDstOutputManager *out = new Fun4AllDstOutputManager("DSTOUT", FullOutFile);
0675     if (Enable::DSTOUT_COMPRESS)
0676     {
0677       ShowerCompress();
0678       DstCompress(out);
0679     }
0680     se->registerOutputManager(out);
0681   }
0682   //-----------------
0683   // Event processing
0684   //-----------------
0685   if (Enable::DISPLAY)
0686   {
0687     DisplayOn();
0688
0689     gROOT->ProcessLine("Fun4AllServer *se = Fun4AllServer::instance();");
0690     gROOT->ProcessLine("PHG4Reco *g4 = (PHG4Reco *) se->getSubsysReco(\"PHG4RECO\");");
0691
0692     cout << "-------------------------------------------------" << endl;
0693     cout << "You are in event display mode. Run one event with" << endl;
0694     cout << "se->run(1)" << endl;
0695     cout << "Run Geant4 command with following examples" << endl;
0696     gROOT->ProcessLine("displaycmd()");
0697
0698     return 0;
0699   }
0700
0701   // if we use a negative number of events we go back to the command line here
0702   if (nEvents < 0)
0703   {
0704     return 0;
0705   }
0706   // if we run the particle generator and use 0 it'll run forever
0707   // for embedding it runs forever if the repeat flag is set
0708   if (nEvents == 0 && !Input::HEPMC && !Input::READHITS && INPUTEMBED::REPEAT)
0709   {
0710     cout << "using 0 for number of events is a bad idea when using particle generators" << endl;
0711     cout << "it will run forever, so I just return without running anything" << endl;
0712     return 0;
0713   }
0714
0715   se->skip(skip);
0716   se->run(nEvents);
0717
0718   //-----
0719   // QA output
0720   //-----
0721
0722   if (Enable::QA) QA_Output(outputroot + "_qa.root");
0723
0724   //-----
0725   // Exit
0726   //-----
0727
0728   se->End();
0729   std::cout << "All done" << std::endl;
0730   delete se;
0731   if (Enable::PRODUCTION)
0732   {
0733     Production_MoveOutput();
0734   }
0735
0736   gSystem->Exit(0);
0737   return 0;
0738 }
0739
0740 # ifndef __CINT__
0741 int main(int argc, char* argv[]) {
0742     if(argc < 1 || argc > 9){
0743         cout << "usage: ./bin/Fun4All_G4_sPHENIX nEvents seed outdir tag outputFile inputFile embed_input_file skip" << endl;
0744         cout << "nEvents: Number of events to generate. Default: 1." << endl;
0745         cout << "seed: random seed to use for reproducibility. Default: 0" << endl;
0746         cout << "outdir: Location of output directory. Default: ." << endl;
0747         cout << "tag: Give a tag to the outputfile. Default: 0" << endl;
0748         cout << "inputfile: Default: https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root" << endl;
0749         cout << "embed_input_file: Default: https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root" << endl;
0750         cout << "skip: Number of events to skip. Default 0." << endl;
0751         return 1;
0752     }
0753
0754      int nEvents = 1;
0755      int seed = 0;
0756      string outdir = ".";
0757      string tag = "0";
0758      string outputFile = "G4sPHENIX.root";
0759      string inputFile = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root";
0760      string embed_input_file = "https://www.phenix.bnl.gov/WWW/publish/phnxbld/sPHENIX/files/sPHENIX_G4Hits_sHijing_9-11fm_00000_00010.root";
0761      int skip = 0;
0762
0763     if(argc >= 2) {
0764         nEvents = atoi(argv[1]);
0765     }
0766     if(argc >= 3) {
0767         seed = atoi(argv[2]);
0768     }
0769     if(argc >= 4) {
0770         outdir = argv[3];
0771     }
0772     if(argc >= 5) {
0773         tag = argv[4];
0774     }
0775     if(argc >= 6) {
0776         outputFile = argv[5];
0777     }
0778     if(argc >= 7) {
0779         inputFile = argv[6];
0780     }
0781     if(argc >= 8) {
0782         embed_input_file = argv[7];
0783     }
0784     if(argc >= 9) {
0785         skip = atoi(argv[8]);
0786     }
0787
0788     Fun4All_G4_sPHENIX(nEvents, seed, outdir, tag, outputFile, inputFile, embed_input_file, skip);
0789
0790     cout << "done" << endl;
0791     return 0;
0792 }
0793 # endif
0794
0795 #endif